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Abstract  of  Dissertation  Presented  to  the  Graduate  Council 
m Partial  Fulfillment  of  the  Requirements  for 
the  Degree  of  Doctor  of  Philosophy 

KINETIC  THEORY  OF  THE  SPATIAL  AND  SPECTRAL 
COUPLING  OF  THE  REACTOR  NEUTRON  FIELD 

By 

Robert  Gifford  Cockrell 
August,  1965 

Chairman:  Dr.  Rafael  B.  Perez 

Major  Department:  Nuclear  Engineering 

A semidirect  variational  method  is  used  to  reduce 
the  integro-dif ferential  form  of  the  general  Boltzmann 
neutron  transport  equation  to  a coupled  set  of  first  order 
ordinary  time-dependent  differential  equations.  The 
independent  variables  in  this  set  of  coupled  equations  are 
the  amplitudes  of  the  neutron  flux  within  each  arbitrarily 
specified  spatial  region  of  the  core  and/or  within  each 
arbitrarily  specified  energy  group  of  the  neutron  spectrum. 

A study  is  made  to  determine  the  accuracy  of  point- 
model  kinetic  theory  in  the  prediction  of  the  fundamental 
prompt  neutron  decay  constant.  The  results  of  the  study 
indicate  that  even  the  1-point-model  kinetics  equation  pre- 
dicts the  decay  constant  to  within  a few  per  cent  of  thd 
exact  value. 
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The  results  of  flux  calculations  for  the  University 
of  Florida  Training  Reactor  are  presented.  The  results 
include  a set  of  18-group  cross-sections  and  graphs  of 
1-  and  2-dimensional  calculations  of  the  flux.  On  the 
basis  of  18-group,  1-dimensional  flux  calculations  the 
effective  delayed  neutron  fraction  is  0.00703  and  the 
mean  prompt  neutron  generation  time  is  247.5  micro- 
seconds. Comparisons  of  theoretical  calculations  with 
experimental  measurements  show  good  agreement. 

In  connection  with  this  work  several  computer  codes 
were  developed  for  use  on  the  University  of  Florida  IBM- 
709  computer.  A brief  discussion  of  these  codes  is  given. 

The  coupled  point-model  representation  of  the 
reactor  neutron  field  should  find  its  greatest  applica- 
biiity  in  the  study  of  transients  and  of  the  propagation 
of  perturbations  in  a reactor  with  weakly  coupled  fuel 
regions. 


CHAPTER  I 


INTRODUCTION 

In  a paper  published  in  the  Proceedings  of  the 
Second  Geneva  Conference  J.  C.  Stewart  made  the  obser- 
vation that  "in  the  analysis  of  certain  types  of  critical 
multiplying  systems,  one  fruitful  approach  to  the  deter- 
mination of  gross  flux  distributions  lies  in  regarding 
the  system  as  made  up  of  several  subcritical  regions 
whose  mutual  interaction  is  sufficient  to  render  the 
system  as  a whole  critical"  (Ste  58) . The  theory 
developed  by  Stewart  was  in  terms  of  the  coupling 
between  regions  containing  fissile  material.  Stewart 
considered  only  the  time-independent  problem,  but  in 
the  same  proceedings  R.  Avery  presented  his  "Theory  of 
Coupled  Reactors"  (Ave  58a)  in  which  he  formulated  a 
method  for  studying  the  kinetics  of  coupled  reactors, 
where  the  term  "coupled"  is  taken  to  mean  that  in  each 
region  of  the  reactor  some  of  the  fission  neutrons  are 
emitted  in  fissions  induced  by  neutrons  born  in  other 
regions . 
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Avery's  theory  presents  a nodal  approach  to  studying 
3pace-dependent  kinetics  problems.  His  method  treats  the 
system  in  terms  of  integral  parameters  which  explicitly 
characterize  the  individual  regions  and  the  coupling 
between  them.  Avery  used  his  theory  of  coupled  reactors 
to  predict  the  kinetic  behavior  of  fast-thermal  coupled 
systems  (Ave  58b) . 

After  studying  Avery's  work,  the  idea  occurred 
to  the  author  that  his  theory  might  also  be  useful  for 
studying  the  kinetic  behavior  of  the  University  of  Florida 
Training  Reactor  (UFTR)  which  is  an  Argonaut-type  reactor 
with  a two-slab  loading  (Len  56,  Dun  58).  Experimentally, 
it  has  been  found  that  a system  of  this  type  exibits  a 
singly  exponential  stable  period  when  off  critical,  but 
that  there  is  a tendency  for  independent  behavior  of  the 
two  slabs  (Bal  59) . When  an  attempt  was  made  to  apply 
Avery's  theory  to  the  study  of  the  kinetics  of  the  UFTR 
several  problems  arose  which  resulted  in  an  abandonment 
of  the  attempt.  Three  of  the  most  serious  problems  are 
summarized  below: 

(1)  The  dependent  variable  in  Avery's 
formulation  is  not  a directly 
measurable  quantity. 

(2)  The  set  of  equations  describing  the 
kinetic  behavior  of  the  coupled 
system  was  formulated  on  the  basis 
of  physical  arguments  and  is  not 
derivable  from  first  principles. 
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(3)  The  system  may  only  be  divided  into 

regions  which  contain  fissile  material. 

A further  search  of  the  literature  revealed  some 
notable  work  by  Dougherty  and  Shen  (Dou  61) . In  this 
work  a semidirect  variational  method  (Hil  52)  is  used 
to  reduce  the  space-  and  time-dependent  neutron  diffusion 
equation  to  a consistent  set  of  ordinary  differential 
equations.  In  this  analysis  the  neutron  flux  is  expanded 
in  a set  of  nonorthogonal  space  modes  with  time-dependent 
coefficients.  This  method  is  more  properly  termed  a 
modal , as  opposed  to  nodal , representation  because  the 
time-dependent  parameters  are  modal  amplitudes  rather 
than  nodal  amplitudes.  There  is  also  a restriction  in 
this  approach.  Each  mode  is  defined  according  to  the  . 
region  where  the  neutrons  represented  by  that  mode  orig- 
inated. Each  mode  extends  over  all  regions  of  the 
reactor.  Since  it  is  the  response  to  local  disturbances, 
either  spatial  or  spectral,  that  is  of  interest  in  coupled 
systems,  a method  is  needed  for  describing  the  coupling 
between  configuration  regions  and  energy  groups  in  addition 
to  the  coupling  between  modes. 

A method  that  has  been  widely  used  (Boy  64;  Sea  64) 

* 

to  study  the  coupling  between  the  fuel  regions  of  weakly 
coupled  cores  is  based  on  Baldwin's  "Kinetics  of  a Reactor 
Composed  of  Two  Loosely  Coupled  Cores"  (Bal  59) . In  this 
method  the  two  regions  of  the  core  are  treated  as  separate 
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subcritical  assemblies  with  a small  interaction.  The 
interaction  between  the  two  regions  is  expressed  very 
simply  by  a term  containing  a factor  e that  gives  the 
fraction  of  the  flux  in  a given  region  at  a time  t that 
was  in  the  other  region  at  a time  t-T , where  t is  the 
delay  time.  The  heuristic  nature  of  this  formulation 
has  served  to  satisfy  the  needs  of  the  experimentalist, 
but  the  ad  hoc  way  in  which  the  parameters  e and  t are 
introduced  places  limitations  on  the  extension  of  this 
work  to  more  complicated  systems  and  makes  e and  t 
difficult  to  predict  theoretically  even  for  this  simple 
system. 

A more  general  nodal  technique  has  been  formulated 
by  Kaplan,  et  al.  (Kap  62) . This  technique  consists  of 
dividing  the  reactor  volume  into  subregions  in  each  of 
which  the  time  behavior  of  some  average  flux  is  deter- 
mined. Nodal  equations  are  obtained  by  integrating  the 
diffusion  equations  over  each  node. 

All  of  the  methods  discussed  above  contain  some 
useful  features,  but  each  contains  undesirable  restric- 
tions. One  notable  restriction  that  has  not  yet  been 
mentioned  is  that  none  of  these  methods  allow  a nodal 
representation  of  the  energy-dependence  as  well  as  the 
space-dependence  of  the  reactor  neutron  field.  In  the 
work  presented  in  this  dissertation  a semidirect  varia- 
tional method  is  used  to  reduce  the  space-,  energy-, 
angle-,  and  time-dependent  integro-dif ferential  form  of 
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the  Boltzmann  neutron  transport  equation  to  a coupled  set 
of  ordinary  first  order  time-dependent  differential  equations. 
This  nodal  representation  of  the  neutron  field  is  so  general 
as  to  include  the  coupling  between  spatial  regions  and/or 
energy  groups,  and  the  flux  used  to  compute  the  coupling 
parameters  may  be  represented  by  any  number  of  eigenmodes. 

All  parameters  appearing  in  this  formulation  are  mathe- 
matically defined  and  are  free  from  ad  hoc  restrictions. 

In  its  simplest  form  this  formulation  yields  the  familiar 
1-point-model  kinetics  equation.  For  a system  repre- 
sented by  one  region,  one  group  and  several  properly 
defined  eigenmodes  (see  Chapter  III) , this  formulation 
is  equivalent  to  that  of  Dougherty  and  Shen  (Dou  61) ; and 
for  a system  represented  by  one  mode,  one  group,  and 
several  regions,  this  formulation  is  equivalent  to  that 
of  Kaplan,  et  al.  (Kap  62). 

To  facilitate  the  utilization  of  the  general 
coupled  kinetics  formulation  presented  in  Chapter  II, 
the  coupling  parameters  are  given  in  terms  of  multigroup 
diffusion  theory  in  Chapter  III.  In  Chapter  IV  an 
analytic  study  of  the  space-dependent  kinetics  of  a 
weakly  coupled  reactor  is  performed  in  order  to  compare 
values  of  the  prompt  neutron  decay  constant  predicted  by 
a space— time  representation,  a coupled  3— point  represen- 
tation, and  a 1-point  representation. 
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Chapter  V presents  a summary  of  the  computer  codes 
that  were  used  to  perform  a computation  of  the  group 
constants , flux,  coupling  parameters,  and  prompt  neutron 
decay  constant  for  the  UFTR.  Chapter  V and  the  Appendix 
contain  tables  and  graphs  of  reactor  physics  parameters 
for  the  UFTR. 


CHAPTER  II 


THEORY 

Introduction 

A semidirect  variational  method  (Hil  52)  is  used  to 
derive  a set  of  first-order  ordinary  differential  equations, 
which  describe  the  kinetics  of  a coupled-core  reactor,  from 
the  integro-differential  form  of  the  general  Boltzmann 
neutron  transport  equation  (Wei  58) . The  independent 
variables  in  this  set  of  coupled  equations  are  the  ampli- 
tudes of  the  neutron  flux  within  each  spatial  region  of 
the  core  and/or  within  each  energy  group  of  the  neutron 
spectrum. 

The  equations  developed  in  this  chapter  should  be 
very  useful  for  studying  the  kinetic  behavior  of  cores 
with  weakly  coupled  regions.  To  aid  in  the  exposition  of 
the  theory  a space-  and  energy-dependent  transfer  function 
is  derived  for  a time  varying  localized  perturbation. 

The  Most  Fundamental  Space-,  Energy-,  and 
Angular-Dependent  Reactor  Kinetics  Equation 

In  terms  of  macroscopic  phenomena  the  fundamental 
expression  describing  the  kinetics  of  the  reactor  neutron 
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field  is  the  Boltzmann  neutron  transport  equation.  This 
equation  describes  the  kinetics  of  the  neutron  field  within 
each  elemental  volume  of  phase  space.  In  terms  of  notation 
similar  to  that  used  by  Weinberg  and  Wigner  (Wei  58)  the 
space-,  angle-,  energy-,  and  time-dependent  integro- 
differential  form  of  the  neutron  transport  equation  is 
given  by: 

1 3 + D 

*(r,E,n,t)  = {x(E)  - l 0 , f . } P*(r,E'  ,$•  ,t)  + 

v (e)  at  1 1 

D 

l Xifi(E)  ci(r,t)  + s(r,E,Q,t)  + 

i=l 

T«Hr,E' ,n*  ,t)  - n»v*(r,E,fi,t)  - 

ET(r,E)  *(r,E,ft,t)  (2.1) 

The  distribution  of  the  precursor  concentration, 
c^(r,t),  for  the  i-th  precursor  group  is  described  by: 

a 

— c.(r,t)  - B.P*(r,E' ,3' ,t)  - X.c. (r,t) 

at-  x 1 1 


i = 1 ,2  , • • • ,D 


(2.2) 
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The  symbols  are  defined  as  follows: 


* (r ,E,n,t) 

= the  angular  flux  (the  number  of  neutrons 
per  unit  area  per  unit  time  per  unit  solid 
angle  per-unit  energy  in  the  volume  element 
dr  abou£  r , traveling  within  the  solid 
angle  dQ  in  the  direction  n,  and  having 
an  energy  in  the  range  dE  containing  E 
at  time  t) 

P*(r,E' ,t) 

“ / dE'  / dfi'  v (E' ) Ef  (r,E' ) ij»(r,E'  ,£$•  ,t) 

e'  n' 

T*(r,E' ,n' ,t) 

■ / dE*  / dn'  E (r,E’-*-E,5'-*-3)  * (r,E*  ,2 ' ,t) 

E-  S* 

(r,E) 

= the  macroscopic  fission  cross-section  for 
neutrons  at  position  r with  energy  E 

v (E) 

= the  mean  number  of  neutrons  resulting 
from  a fission  caused  by  a neutron 
with  the  energy  E 

Eglr^'+E^'-*^) 

= the  macroscopic  scattering  cross-section 
at  position  r^for  changing  the  energy  and 
direction  E' ,n'  into  an  energy  and 
direction  range  dE,  dn  about  E,n 

X(E) 

“ the  fission  spectrum  for  all  neutrons, 
both  prompt  and  delayed 

^(E) 

= the  energy  frequency  function,  or 
spectrum,  of  the  neutrons  emitted  by 
the  i-th  precursor  group 

ST(r,E) 

= the  total  macroscopic  collision  c|oss- 
section  for  neutrons  at  position  r 
with  energy  E 

(r,t) 

= the  concentration  of  the  i-th  precursor 
group  at  point  r and  time  t 
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* the  decay  constant  of  the  i-th 
precursor  group 


0^  = the  fraction  of  fission  neutrons 

that  are  emitted  by  the  i-th  pre- 
cursor group 


^•7  = th<g  leakage  operator  (operates  on 

i(i(r,E,0,t)  to  give  the  rate  of  change 
of  the  angular  |lux  with  respect  to 
space  variable  r in  the  direction  n) 


s(r,E,Q,t)  = the  number  of  neutrons  (per  unit  time 
per  unit  energy  per-unit  solid  angle 
per  unit  volume  at  r,E,  and  n)  at  time 
t from  an  external  source 


Implicit  in  the  formulation  of  Eqs.  (2.1)  and  (2.2) 
is  the  assumption  that  the  precursors  do  not  diffuse, 
which  limits  this  treatment  to  fixed  fuel  reactors.  It 
appears  that  Eqs.  (2..1)  and  (2.2)  are  also  limited  to 
reactors  with  only  one  fissile  isotope;  however,  reactors 
with  several  fissile  isotopes  are  included  if  the  production 
operator  is  written  in  the  forms 

x(E)P*(r,E',fi'  ,t)  = l Xj(E)  / dE'  / dQ*  v ^ (E ' ) • 

j e*  if' 

, < 

E^(r,E' ) 4>  (r,E'  ,t) 

and 

eifi(E)P*(r,El  ,t)  = l eJfJ(E)  / dE'  J dO'v^E')  • 

j E'  ft' 

I^(r,E') * ( r ,E ' ,ft' ,t) 
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where  j refers  to  some  j-th  fissile  isotope,  and  the  i-th 
precursor  group  is  assumed  to  be  the  same  for  all  species. 

The  boundary  conditions  for  Eqs.  (21.)  and  (2.2)  are 
the  usual  ones  (Meg  60)  but  are  stated  here  for  later 
reference : 

(1)  The  flux,  and  likewise  the  precursor  concen- 
trations, must  be  finite  and  nonnegative  for 
all  points  in  phase  space  except  possibly  at 
a finite  number  of  source  points 

(2)  The  flux  ijj  (r+efj  ,E , 0 , t)  must  be  a^continuous 
function  of  e for  all  values  of  0,  E,  and  t 

(3)  At  a vacuum  interface  (the  outermost  boundary 
of  a system)  there  must  be  no  neutrons  return- 
ing to  the  reactor  from  any  direction  in  the 
vacuum 

Because  of  the  complexity  of  Eqs.  (2.1)  and  (2.2) 
the  problems  for  which  exact  solutions  may  be  found  are 
few  in  number  and  apply  to  ideal  systems.  In  order  to 
analyze  problems  that  arise  in  practice  it  is  necessary 
to  make  simplifying  assumptions.  In  the  study  of  reactor 
kinetics  one  such  assumption  is  that  the  neutron  field  has 
the  same  time-dependence  in  each  elemental  volume  of  phase 
space  for  all  values  of  r,E,Q.  This  assumption  permits 
one  to  write  the  function  representing  the  angular  flux 
in  a separated  form  as  follows: 

♦ (r,E,$)T(t) 


*(r,E,(5,t) 


SB 


(2.3) 
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There  are  many  methods , both  numerical  and  analytical  for 
obtaining  approximate  solutions  to  the  steady  state  flux, 
<J>(r,E,n).  The  steady  state  flux  along  with  its  adjoint 
$ (r,E ,ft)  is  then  used  to  reduce  the  transport  equation  to 
a first  order  ordinary  differential  equation  with  T(t)  as 
the  dependent  variable.  This  equation  is  the  familiar 
point-model  reactor  kinetics  equation  (Shu  61)  . The 
derivation  of  the  point-model  equation  is  discussed  in 
great  detail  by  Henry  (Hen  55)  and  therefore  will  not  be 
given  here.  In  the  derivation  of  the  coupled  kinetics 
equations  given  in  the  next  section  the  assumption  of 
separability  is  less  restrictive  since  it  only  requires 
the  separability  of  time  and  phase  space  over  a subspace 
of  the  domain  containing  the  reactor  neutron  field. 

Derivation  of  the  Coupled  Kinetics  Equations 

The  derivation  which  follows  is  based  on  a semi- 
direct  variational  method  (Hil  52) . This  method  consists 
of  varying  the  functional  after  substituting  and  inte- 
grating the  trial  functions  over  the  domain  of  all  but  a 
single  variable.  It  is  assumed  that  the  distribution  of 
the  angular  flux  is  known  at  steady  state.  The  purpose  of 
using  the  variational  method  is  to  determine  the  "best" 
approximating  function  for  the  time-dependence  of  the  flux. 
By  definition  the  "best"  approximating  function  is  that 
trial  function  which  makes  the  functional  a minimum.  A 
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problem  arises  regarding  the  minimizing  of  the  functional 
because  the  variational  principle  is  not  a maximum  or 
minimum  principle  but  rather  it  is  an  extremum  principle. 
Sometimes  it  is  possible  to  study  the  second-order  varia- 
tion of  the  functional  to  determine  if  the  trial  function 
produces  a maximum  or  a minimum,  but  this  approach  is 
extremely  difficult  for  very  complicated  problems.  This 
fact,  however,  is  of  no  importance  as  regards  the  deri- 
vation of  the  equations  of  motion,  since  only  the  extremum 
condition  is  used  (Lan  60).  Usually,  it  is  possible  to 
tell  from  the  physical  situation  whether  the  extremum  is 
a maximum  or  a minimum. 

In  the  derivation  which  follows  it  should  be  noted 
that  the  presentation  will  be  entirely  formal;  any  questions 
of  existence  or  uniqueness  are  assumed  to  be  resolved  on 
the  basis  of  physical  considerations.  Selengut's  "Varia- 
tional Analysis  of  Multi-Dimensional  Systems"  (Sel  59)  is 
an  excellent  reference  for  a more  detailed  discussion  of 
the  variational  principle  as  applied  to  reactor  physics. 

Equations  (2.1)  and  (2.2)  may  be  written  in  the 
form  (Per  62) : 

i 3 

h?  - v-i  — y - s 

3 1 


(2.4) 
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where 

m 

f =*  [i|»  (r,E,fl,t)  c1(r,t)  •••  c..(r,t)  •••  cD(r,t)]  (2  5) 

«* 

I * [s(r#E,S,t)  0 •••  0 •••  0 ] (2.6) 


The  symbol  indicates  the  transpose  of  the  vector. 


H = 


D 


{X(E>-  l 6ifi } | + T-fi.V-Z, 
i=l 


81f1P/X 


A1F1 


0 


• • • u • • • 


“XDfD 


(2.7) 

The  parameter  A has  been  inserted  in  the  matrix 
represented  by  H to  insure  the  existence  of  a solution  for 
the  homogeneous  case.  In  general  there  is  an  infinite  set 
of  the  eigenvalues  represented  by  A with  a corresponding 
set  of  eigenfunctions. 

The  symbol  V ^ represents  the  diagonal  matrix 
given  by 


_1  A 
V * diag  [ — f 


, • • • f , • • • f ] 
1 i DJ 


(2.8) 


The  Lagrangian  formulation  of  a variational 
principle  for  the  equation  given  by  Eq.  (2.4)  is: 
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L(*V) 


/ dt  / dr  / dE  / dfrr+ 
° R E 0 


3 

— Y-HY-3] 
3t 


(2.9) 


where 


*►  "f* 

» 


C'F+  (r,E,0,t)  c*(r,E,t) 


•••  c|(r,E,t) 


c*(rfE,t) 


The  condition  that  the  first  variation  of  the 
functional  ,4')  be  an  extremum  with  respect  to  small 

variations  in  the  functions  and  I yields  an  expression 
for  the  weighting  function  as  follows: 

6L(6Y+,SY)  » / dt  / dr  / dE  / dft  { 6$+  [V-1  — + 

° R E ft  9t 

V 

+ , a 

t CV"1  — fif  - H6f]} 

(2.10) 


By  definition  of  the  adjoint  operator  (Sel  59)  the 
last  term  in  Eq.  (2.10)  may  be  rearranged  to  yield 

(?+,H6f)  = (5f,H+?+)  (2 

1* 

where  H is  the  adjoint  of  the  matrix  H.  The  parenthesis 
notation  indicates  integration  over  phase  space.  The 
adjoint  operator  H+  is  given  by: 
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D 


[ {X (E) - l 0.f.}  |]+  + 
i=»l 


T + +fi*  7-Zr 


H' 


Xlfl 


XDfD 


C61f1P/X]+  ...  '[8DfDP/X]  + 


'Xlfl 


o ...  -X  f 

D D 


(2.12) 

where 

f ° p + + v(E)Ef(r,E) 

C { X ( E ) — l e.f.(E)}  / = / dE*  / d3'  . 

i=l  x 


D 

{ X (E ' ) “ l 0ifi(E')  }*+(r,E'  ,3'  ,t) 

i=l  (2.13) 


and 


T%+(r,E'  ,n*  ,t)  » / dE'  / dfl'  E (r,E-E'  ,M')*  + (r,E'  ,3*  ,t) 

E'  3' 


17 


An  integration  by  parts  of  the  next  to  last  term 
of  Eq.  (2.10)  yields: 


00  3 

/ dt  / dr  / dE  / dfi  $+  v”1  — 

0 R E 3 3t 


/ dr  / dE  / ds$  • 
R E S 


[1'  + V”15l]" 

t=  0 


(6?, 


(2.14) 


The  requirement  that  the  first  variation  of  f vanish  at 
the  end  points  of  the  time  domain  causes  the  first  term  on 
the  right  hand  side  of  Eq.  (2.14)  to  vanish.  Inserting 
Eqs.  (2.11)  and  (2.14)  into  Eq.  (2.10)  yields): 

4L(5^  ,5$)  = / dt  / dr  / dE  / dft  [V-1  

0 R E 0 3t 


5 ¥ [V-1  — $+  + H+?+]} 

3 1 

The  condition  for  an  extremum  is  that 


(2.15) 


6L(6yf,6f)  » 0 

which  requires  that  the  integrand  of  Eq.  (2.15)  be  zero. 
Since  the  variations  6Y  + and  6$  are  arbitrary  and  in 
general  are  not  equal  to  zero  except  at  the  boundaries 
of  the  phase  space  spanning  the  reactor  neutron  field, 
then  the  coefficients  of  6Y  + and  6?  must  vanish.  This 
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requirement  yields  Eq.  (2.14)  and  Eq.  (2.16). 

(2.4)  is  repeated  here  for  convenience. 

. , 3 

E'V  = V I - § 

3 1 

5i+ , - v-1  -1 

3 1 


Equation 


(2.4) 


(2.16) 


Thus,  any  trial  function  used  to  approximate  ¥ will  make 

the  value  of  L(Y+,$)  an  extremum  if  the  weighting  function 

^ satisfies  Eq.  (2.16)  and  the  same  boundary  conditions  as 

This  latter  point  raises  an  interesting  question  re- 

garding  the  boundary  conditions  to  be  placed  on  the  time 

variable.  The  first  order  time  derivative  in  Eqs.  (2.4) 

and  (2.16)  permits  only  one  boundary  condition  with  respect 

to  the  time  variable  for  each  component  of  the  vectors  $ 

and  4'+.  However,  the  requirement  that  $+v“1Sy  vanish  at 

both  limits  of  the  time  domain  seems  to  imply  that  I is 

known  at  t=0  and  also  at  t=»#  it  will  be  shown  that  this 

is  not  correct,  but  that  the  condition  that  $+v"16$|*° 

t=*0 

vanish  is  not  valid. 

, «* 

Multiplying  Eqs.  (2.4)  by  V and  Eq.  (2.16)  by  V, 

integrating  over  all  values  of  t,E,ft  and  r,  and  subtracting 
Eq.  (2.16)  from  Eq.  (2.4)  yields: 


(y+,H?)  - (?,Htft)  = (f+,v“1  — Y)  + (f,V_1  -1  f+)  - (*+,3) 

at  3 1 


a - 


(2.17) 


zero ; 
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The  left-hand  side  of  Eq.  (2.17)  is  identically 
rearranging  the  remainder  of  Eg.  (2.17)  yields: 


— (f+V_1Y)  = (|+,£) 

at 


(2.18) 


Eguation  (2.18)  indicates  that,  if  there  is  no  external  • 

source,  the  product  is  invariant  with  respect  to  time. 

In  other  words,  the  adjoint  weighted  flux  is  conserved.' 

Therefore,  it  is  not  necessary  to  know  the  flux  precisely 

at  both  time  limits  to  insure  that  ^+v”1jl|”  vanishes. 

t=0 

If  a source  is  present  during  some  specified  period  of 
time,  then  the  change  in  the  adjoint  weighted  flux  is 
given  by  Eg.  (2.18) . 

There  are  numerous  methods  for  obtaining  time- 
independent  approximations  of  the  solutions  for  the  source- 
less, steady-state  form  of  Egs.  (2.4)  and  (2.16).  Therefore, 
in  the  derivation  which  follows,  a time-independent  solution 
is  assumed  to  be  known.  The  trial  function  to  be  used  in 
the  functional  given  by  Eg.  (2.9)  is: 

M M 

* " l *n(r,E,$,t)  * l '?n(?,E,fi)N]Jn)  (t) 

n=1  n=1  Y (2.19) 

and 

M . M 

- l = l Tjl(?,E,5)S|Jra>t(t) 

m*l  m=l 


Y 


(2.20) 
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where 


N^n) (t) 


the  amplitude  of  $n(r,E,n,t)  in 

region  and  energy  range  for 

any  direction  the  configuration 

space  has  been  divided  into  K regions 

R1 'R2 ' * * * 'Rk ' * * * ,RK  and  the  ener?y 
space  has  been  divided  into  r groups 


E1'E2'*‘* 


#E. 


rE, 


The  components  of  N^n  (t)  and  (t)  are  given  by: 

Y Y 


'^n)  (ti 

Y 


l^n)  (t)  c^(t) 

Y ' Y 


Y 


(2.21) 


N<m)  (t) 


[N^ml  (t)  (t) 

Y 'Y 


c«  itn 

Y 


(2.22) 


and  the  initial  conditions  are: 

N,[n)  (0)  = [ 1 1 •••  1 ] 

Y 

f 

Nj[m)  (0)  = [ 1 1 •••  1 ] 

Y 

for  all  values  of  k and  y.  The  quantities  Yn(r,E,fi)  and 
2)  are  diagonal  matrices  given  as  follows: 
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= diag  [<J>  (r,E,S$)  c.  (r)  c~  (r)  •••  c (r)  ] 

11  n J.  fin  zfn  D ,n 

(2.23) 

,E ,«)  - diag  [**(r,E,fi)  c^m(r,E)  c^m(r,E)  ... 


CD,m(r'E)] 


(2.24) 


The  quantities  4'n(r,E/n)  and  T + (r,E,j$)  satisfy  the 
following  equations  and  the  boundary  conditions  specified 
earlier  j 


H ¥ 5 = 0 

n n 


H+4'+t  = 0 

m m 


where 

t - [ 1 1 •••  1] 


(2.25) 

(2.26) 


and  is  given  by  Eq.  (2.7)  for  X=X^  and  is  given 
by  Eq.  (2.12)  for  X=X ^ . The  quantity  X ^ is  the 
eigenvalue  corresponding  to  the  i-th  eigenfunction  of 
Eqs.  (2.25)  and  (2.26). 

The  symbols  used  in  Eqs.  (2.21) -(2.24)  have  the 
following  definitions: 

■>  *f. 

<j>n(r,E,Q)  = the  n-th  eigenmode  of  the  angular 

flux  at  steady  state 

<f>^(r,E,^)  = the  adjoint  of  <j>n(r,E,n) 
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°i,n(?) 


ci,n(r'E) 


N^n) (t) 


N^n)  (t) 


C!k‘  <« 

Y 


ciS}  (t) 


= the  n-th  eigenmode  of  the  precursor 
concentration  distribution  for  the 
i-th  precursor  group 

= the  adjoint  of  ci  n(r).  (Referring 

back  to  Eqs . (2.7)  and  (2.12),  it  is 

observed  that  ct  n is  energy-dependent 

although  c.  is  not) 
l ,n 

= the  amplitude  of  the  n-th  eigenmode 
of  the  angular  flux  in  region  R, 
and  group  at  time  t K 

= the  amplitude  of  the  n-th  eigenmode 
of  the  adjoint  angular  flux  in  region 
R^  and  group  E^  at  time  t 

= the  amplitude  of  the  n-th  eigenmode 
of  the  i-th  precursor  group  concen- 
tration in  region  R.  and  group  E 
at  time  t K Y 

= the  amplitude  of  the  n-th  eigenmode 
of  the  adjoint  to  the  i-th  precursor 
group  concentration  in  region  R.  and  • 
group  E at  time  t K 


Before  substituting  the  trial  function  into  Eq.  (2.9), 
the  functional  is  written  in  the  form: 


L($+,¥) 


/ dt  l l / dr  / dE  / dn  4'  + [v“1  — 

0 k=l  y=l  R.  E $ 

k y 

(2.27) 


In  order  to  simplify  the  mathematical  expressions  the 
following  parenthesis  notation  will  be  used: 
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(Gg,Ff)k  = f dr  f dE  f dn  [Gg]  [Ff] 

Y *k  Ey  5 ( 

where  G is  operating  on  g and  F is  operating  on  f.  If 
either  F or  G is  not  specified,  then  it  is  taken  to  be 
the  identity  operator. 

Substituting  the  trial  function  and  its  adjoint 
given  by  Eqs.  (2.19)  and  (2.20)  into  Eq.  (2.27)  yields: 


„ K T M N 3 

MT\*>  - / at  £ l l l (F.V1  -ln)k 

0 k=l  Y=1  m=l  n=l  3t  Y 


(Y  , HY  ), 
m'  n k 
Y 


( u? ^ <?[) 

1 Vb,k 

Y 


(2.29) 


Assuming  separability,  within  region  R,  and  group  E , 

K Y 

yields  for  the  first  term  in  the  integrand: 


(fJ.V-1  — Y ) 
m 3 1 n y 


S£m)  (t)  ^m,n)  i_^n)(t) 

Y Y dt  Y 


(2.30) 


where  N^11^  and  N^1"^  are  given  by  Eqs.  (2.21)  and 
Y Y 

(2.22)  and  is  a diagonal  matrix  given  by 

Y 


A(m,n)  ^ [(<)i+,  i $ ) (ct  ),  ••• 

K ^ rra'  v n k l,m'  1 l,n'k 

Y Y Y 


(cD,m' fDCD,n) k] 
Y 


(2.31) 
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In  a similar  manner  the  third  term  in  the  integrand 
of  Eq.  (2.29)  may  be  written  in  the  form: 


(y+,$) 


k 

Y 


Y 


where 

Y 


is  a vector  with  (D+l) 


components : 


(2.32) 


[ (♦m'S)k  ^ ’**  ® ***  0 ] 

Y 


(2.33) 


The  symbol  "s"  in  Eq.  (2.33)  was  defined  earlier  and 
represents  an  external  source  of  neutrons. 

The  second  term  in  the  integrand  of  Eq.  (2.29)  is 
somewhat  more  complicated  than  the  other  two  terms.  It 
is  advantageous  to  consider  it  term  by  term.  The  only 
non-zero  elements  in  the  matrix  H are  the  following  ones: 


h 


1,1 


D 

* T - I «ifi  f + T - 

i=l 


(2.34) 


h.  ...  = X . f . 

1,  l+l  l i 


h. 


i+i,  i ■ 


hi+l,  i+X  " _xifi 


for  i = 1 ,2 , • • • ,D. 


(2.35) 

(2.36) 

(2.37) 


• • • 
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The  only  non-zero  elements  of  the  matrix  represented 
by  are  given  by 


element  (1,1)  = 


(xt+,  ^ * ), 
m'  X n k 
y 


l VV"'  f *n>k 
i-1  * 


<V  T?n>k  - <V  5'vJn>k  ’ 'VVn'k 

y y y 


(2.38) 


element  (1,  i+1)  = X. (cT  ,f,c.  ). 

i l ,m  l i ,n  k 
Y 


(2.39) 


element  (i+1,  1)  = 6^**,^  ~ 5)k 


(2.40) 


element  (i+1,  i+1)  =* 


-X  . c?  ,f  .c.  ). 

l i ,m'  l l ,n  k 
' Y 


(2.41) 


for  i=l , 2 , • • • ,D. 

A bar  over  a function  in  Eqs.  (2 . 38) -(2 . 41) 
indicates  that  the  function  is  still  time-dependent;  that 
is,  separability  has  not  yet  been  assumed.  Consider  now 
the  first  term  of  Eq.  (2.38): 

(xJm'  r Vk  = / / dE  / dS$  x(E)*^(r,E,^,t)  / dE 1 • 

R,  E 3 E ' 

k Y 

/ dn  v(E') Ef (r,E')*n(r,E' ,t)  = 

Q 
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l (X*;.  f %>k  N^>(t) 

«=1  y X y a (2.42) 


where 


pQ$n  ■ / dE'  / dfi ' v(E')If (r#E')*n(r,E^') 

E*  i?' 

a 


(2.43) 


It  will  prove  useful  later  to  write  Eq.  (2.42)  in  the 
following  form: 


(**m'  T Vk 

y 


G 

i 

a=l 


N 


(m) 

k 

Y 


(*m'  I Vk  , , 

(t)  1 N*n)  (t) 

. (m,n)  a 

k 

a+y 


(2.44) 


where 


. (m,n) 
vk 
a->  y 


<•!'  b Vk 


(*‘m'  ~ Vk 

x y 


(2.45) 


The  parameter  A^m'n^  refers  to  neutrons  existing  in  the 

a+y 

n-th  eigenstate  in  region  and  group  Eq.  It  gives  the 
average  time  for  transfer  from  group  to  group  E^  per 
neutron  in  group  by  means  of  regeneration.  The 
averaging  process  is  weighted  with  respect  to  the 
importance  of  this  transfer  to  the  m-th  eigenstate. 
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The  second  term  in  Eq.  (2.38)  is  similar  to  the 
first  and  may  be  written  in  the  form: 


6i(Vl'  f Vk 

y 


G 

l 

a=l 


Nj[n)  (t) 


(m,n) 

ei,k 

a-»V 

, (m,n) 

Ak 

a->y 


<*m'  k Vk 

Y 


^n)  (ti 


(2.46) 


where 


o (m»n) 

i*y 


(f.$  ,P  $ ),  8. 

i m'  a n k l 


(X$+,P  $ K 

' A m'  a n k 


(2.47) 


and 


8 


(m,n) 

k 

a-*y 


D 


i=l  «-*-Y 


(2.48) 


The  parameter  8^m'n^  gives  the  effective  delayed  neutron 

"a+y 

fraction  averaged  over  the  spectrum  of  the  neutrons  in  the 
n-th  eigenstate  in  region  that  are  transferred  from 
group  E^  to  group  E^  by  regeneration. 

The  third  term  in  Eq.  (2.38)  follows  in  similar 
fashion : 


($+,T®  ). 

' m'  n k 

Y 


r 

I 

a = l 


N 


(m) 

k 

Y 


(m) 

a n'k  “k 
Y a 


(4>+,T  $ ),  N 
m'  ~ " 


(2.49) 


28 


r 

I 


N 


(m) 

k 

Y 


<‘m'  h Vk 


(ro,n) 


k 

a-*  y 


N 


(n) 

k 

a 


where 


_ (m,n) 
Tk 
a-^y 


(tm-  7 Vk 


(*  ,T  $ ) , 
m'  a n k 


and 


(2.50) 


(2.51) 


T ♦ =•  / dE’  / d$'  E (r,E'-»-E,sWfi)  • (r,E'  ) 

a n J ' s n 

E'  O' 

a 


(2.52) 


The  parameter  T^m'n^  refers  to  neutrons  existing  in  the 

a+Y 

n-th  eigenstate  in  region  and  group  Ea.  It  gives  the 
average  time  for  transfer  from  group  E^  to  group  E^  per 
neutron  in  group  E^  by  means  of  scattering  collisions. 

The  fourth  term  in  Eq.  (2.38)  represents  the 
following  expressions 


($+,fi*V$  ) = / dr  / dE  / do  (r  ,E  , $ , t)  (r  ,E  ,2  , t) 

m n k . ^ ** 

Y R,.  E 3 


(2.53) 
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Applying  the  divergence  theorem  (Spi  59)  to  the  evaluation 
of  the  integral  over  the  space  variable  yields 

($*,n«V$n)k  = J dE  J dQ  J dsknk*ft  $^(sk,E,ft,t)  $n(sk,E,3,t)  - 

Y ^ n sk 

/ dr  / dE  / dfi  **(r,E,fi,t)fi»v**(r,E,fi,t) 

Rk  Ey  (2.54) 

where  Sk  is  the  surface  enclosing  region  Rk,  sk,  is  the 

value  of  r on  the  surface  Sk,  and  nk  is  the  outward  drawn 

unit  normal  to  region  Rk>  The  vector  (sk,E  ,j^,t)  gives 

the  neutron  current  in  the  direction  n.  If  n is 

9 

positive,  the  neutrons  are  leaving  region  Rk;  if  nk*S$ 
is  negative,  the  neutrons  are  entering  region  Rk.  With 
this  consideration  in  mind,  the  first  term  of  Eq.  (2.54)  may 
be  expressed  as  the  sum  of  a loss  term  and  a source  term 
as  follows: 


/ dE  / d2  / dsk  (sk,E,S,t) 


ey  2 sk 


J dE  / d2  J dsk  nk-3*^(sk,E,S,t)*n(sk,E,2,t)  - 

2 sk 


K 


l I dE  / d^  / dsk  nk*^«^(sk,E,^ ,t)  ?n(sk,E,j^ ,t) 
j = l E S 


y j k,  j 


(2.55) 
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where 


that  part  of  the  surface 
common  interface  with  the 


Sk  that  has 
surface 


those  values^of  0 for  which  the 
dot  product  n..*fi  is  positive 


a 


Inserting  Eq.  (2.55)  into  Eq.  (2.54)  and  rearranging 
terms  yields  the  following  expression  for  the  last  two 
terms  of  Eq.  (2.38): 
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where 
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l N 

j“l 

j 


(m) 


( <t>^ , i $ ) 
tnf  v n'k 


(t) 


(m,n) 

j-k 

Y 


— N!n)  (t) 
3 


(2.56) 


(m,n)  _ 
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(2.57) 
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( . i.  $ ) 

m'  v n'k 


/ dE  J / ds,  n • $<&*$ 

1 ' k k m n 


E. 


3k , j 


(m,n) 


(2.58) 


The  quantity  gives  the  average  inverse  loss 

Y 

rate  of  neutrons  belonging  to  the  n-th  mode  in  region 

and  group  E^.  The  quantity  gives  the  average  inverse 

Y 

transfer  rate  from  region  R_.  to  region  Rj^  for  neutrons 
belonging  to  the  n-th  mode  in  group  E^ . 

Collecting  the  terms  expressed  by  Eqs.  (2.44),  (2.46), 
(2.50)  and  (2.56)  and  substituting  them  into  Eq.  (2.38) 
yields  the  following  expression  for  element  (1,1)  of  the 

(?m'"fn>k  matrlx! 

Y 
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element  (1,1) 
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(2.59) 


Element  (1,  i+1)  expressed  by  Eq.  (2.39)  may  be 
written  in  separated  form  directly: 


element  (1,  i+1)  * cfra^  (c.  ,f,Cl  ) x C 

' i,k  v i ,m' ri ci ,n;  k i ci,k 


(m) 

Ai" 

Y Y Y 

for  i “ 1,2 , • • • ,D. 

Element  (i+1,  1)  expressed  by  Eq.  (2.40)  must  be 
manipulated  as  follows: 


(2.60) 
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element  (1,  i+1)  =»  N 


(m) 
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r Bi,k 
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I 

a=l  . (m,n) 
k 

a+y 
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m'  v n k k 
Y a 


(2.61) 


for  i = 1 ,2 , • • • ,D. 

The  diagonal  terms,  element  (i+1,  i+l)f  expressed 
by  Eq.  (2.41)  are  given  by: 


element  (i+1,  i+1) 


for  i = 1,2,  •• • , D. 


\ c(ra)  (c^  f c ) p(n) 
i i,k  ci ,m' £ici ,n  k Ci,k 

Y Y 


(2.62) 


. Collecting  the  terms  given  by  Eqs.  (2.30),  (2.32), 
(2.59),  (2.60),  (2.61)  and  (2.62),  rearranging  and  inserting 
them  into  the  expression  for  the  Lagrangian  given  by 
Eq.  (2.29)  yields: 


L(N+,N,N) 
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(m) 
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K 


I D*"'"1  A^"'n>  5<n>(t>  + L<ra'">  ^n,(t)( 

j=l  y y y y y y 


(2.63) 


where 


N1 

-► 

N 

3 


represents  the  vectors  N*  ^ (t) 

k 

Y 

represents  the  vectors  N^n) (t)  and 

d j*  Y 

= 3tN 


The  matrix  A^m'n)  is  given  by  Eq.  (2.29).  The 


other  matrices  are  defined  as  follows: 


B 


(m,n) 

k 

a+y 


a-+y 


. (m,n) 
Ak 
a+y 


fl(m,n)  . 
Bl,k  ' 
a-*-y 


BD,k  ' 
a-*y 


t (m,n) 
k 

a+y 


, (m,n) 
vk 
a-*-y 


(m,n) 

vk 

*-y 


(2.64) 
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n(m,n) 
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(2.65) 
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The  necessary  condition  that  the  N+  and  N make 

the  Lagrangian  an  extremum  is  that  the  Euler-Lagrange 

equations  (Mor  53)  be  satisfied,  i.e., 

% 


• 0 • 

F(ft+  /$/&)  - — [ F(ft+ ,$,$)] 


dt 


3N 


- 0 


(2.67) 


and 


— F(i$+,^,N) 
3N 


d 

dt 


— 7 F(N+  ,N,N)  ] - 0 

3N 


(2.68) 


where  F(N+,N,N)  represents  the  integrand  of  the  functional 
expressed  by  Eg.  (2.63) . 

Application  of  Eq.  (2.67)  yields  the  equations  of 
motion  for  the  components  of 
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j“l  Y Y 
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(2.69) 

1 . 2 ,  • • • , M mode  s 

1.2 , • • • ,r  energy  groups 

1 .2 ,  • • • ,K  regions 
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The  equations  of  motion  for  the  components  of 
may  be  found  with  Eq.  (2.68) ; however,  since  these  equations 
are  not  of  much  practical  use,  they  will  not  be  given  here. 
Let  us  proceed  with  an  examination  of  Eq.  (2.69). 

Equation  (2.69)  will  be  clarified  somewhat  if  the 
set  of  equations  represented  by  the  matrix  notation  is 
written  out  in  full  as  follows: 
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1 , 2 , • • • ,M  modes 

1 .2 ,  • • • ,K  regions 

1.2, « • • ,r  energy  groups 

1 . 2 , • • • ,D  precursor  groups 
the  Kronecker  delta 


When  the  kinetics  of  the  neutron  field  is  represented 
by  one  mode,  one  region,  and  one  energy  group,  Eqs.  (2.70) 
and  (2.71)  reduce  to 


d (1-6)  N ( t)  D (cT,fici) 
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(2.73) 
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where 

A = generation  time 

l = (1/6  - 1/t)  ^ = lifetime 

and  if  the  normalization  factor  (c? , f ici) / ( $+ , i $)  is 
absorbed  into  (^(t),  then  it  is  seen  that  Eqs.  (2.72)  and 
(2.73)  are  the  familiar  point-model  reactor  kinetics 
equations.  This  result  is,  of  course,  as  it  should  be; 
however,  in  connection  with  this,  an  important  observation 
should  be  made.  Note  that  the  method  used  to  derive  the 
very  complicated  coupled  kinetics  equations  also  yields 
the  point-model  kinetics  equation.  Now,  the  simplicity  of 
the  point-model  kinetics  equations  are  appealing  in  com- 
parison with  the  appallingly  complicated  coupled  kinetics 
equations  expressed  by  Eqs.  (2.70)  and  (2.71).  However, 
the  idea  is  not  to  complicate  a simple  problem  but  rather 
to  simplify  a complicated  kinetics  problem  that  is  not 
well  represented  by  Eqs.  (2.72)  and  (2.73).  Equations 
(2.70)  and  (2.71)  should  find  their  greatest  applicability 
in  the  study  of  reactor  cores  with  weakly  coupled  fuel 
regions . 

The  derivation  of  the  coupled  kinetics  equations 
would  not  be  complete  without  showing  how  to  reconstruct 
the  flux  in  terms  of  the  coupled  parameters.  The  components 
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of  the  vector  4*  were  written  in  separable  form  according 
to  Eq.  (2.19)  : 

M 

f , I *n(r,E,ft)  Nj|n)  (t) 

n=1  y (2.19) 

where  V is  defined  by  Eq.  (2.5);  *n(r,E,fi)  is  a matrix 
of  the  steady  state  values  of  the  angular  flux  and  the 
precursor  concentrations  and  is  defined  by  Eq.  (2.23); 
an(3  (t)  is  the  vector  defined  by  Eq.  (2.21).  Multi- 

Y 

plying  Y by  'F^(r,E,fi)V  \ integrating  over  region  R^, 
all  0,  and  energy  group  E^ , and  summing  over  all  values 
of  m yields: 
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I I 
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tT«'v'lT„»k  skn)(t) 


(2.74) 


Assuming  the  eigenfunctions  are  members  of  a 
complete  set  and  that  M such  modes  are  adequate  to  describe 
the  functions  represented  by  y*,  the  left  hand  side  of 
Eq.  (2.74)  is  given  by: 

M 

I = A* 

m=l  Y Y Y (2.75) 

where  A^  is  the  matrix  defined  by  Eq.  (2.31)  summed  over  all 
Y 

values  of  m and  n. 
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On  the  left  hand  side  of  Eq.  (2.74) 

(,'m'V~1Vk  ’ ^m'n> 

y y (2.76) 

Substituting  Eq.  (2.75)  and  (2.76)  into  Eq.  (2.74)  and 

operating  from  the  left  with  a"1  yields: 

Y 

M M 

\(t)  - I { l A"1  A,^'n)  } t^n)  (t) 

Y n=l  m**l  Y y Y (2.77) 

Equation  (2.77)  shows  that  the  amplitude  of  the  neutron 
population  in  region  and  group  is  given  by  a weighted 
sum  of  the  amplitudes  of  the  neutron  population  in  each  mode. 
The  term  in  {}  of  Eq.  (2.77)  gives  the  fractional  contri- 
bution of  the  n-th  mode  to  the  total  population  in  region 
and  energy  group  E^ . 

Similarly , it  can  be  shown  that  the  time-dependence 
of  the  flux  average  over  the  entire  reactor  is  given  by 

K r M M 

ff<t)  =1111  A-1A|Jm'n) 

k=l  n=l  m=l  y 

Uhere  M N 

A - I I 

m=*l  n=l 

To  aid  in  the  exposition  of  the  theory  a method 
for  determining  a space-  and  energy-dependent  transfer 


^n)  (ti 

Y (2.78) 
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function  is  discussed  in  the  next  section.  To  be  specific 
the  pile  oscillator  experiment  is  discussed;  however,  the 
resulting  expressions  are  also  applicable  to  other  studies 
involving  a localized  time-dependent  perturbation  in  the 
absorption  cross-section. 

Application  to  the  Pile  Oscillator  Experiment 

One  method  of  determining  the  zero-power  transfer 

function  for  a reactor  is  by  means  of  a pile  oscillator 

experiment.  In  such  an  experiment,  a time-dependent 

perturbation  is  created  at  some  point  in  the  reactor 

while  fluctuations  in  the  neutron  population  are  measured 

with  a detector  at  some  other  point.  The  ratio  of  the 

% 

detector  output  to  the  input  perturbation  yields  the 
transfer  function.  In  a reactor  with  a compact,  or 
"tightly  coupled"  core,  such  an  experiment  yields  meaning- 
ful results;  however,  in  composite  or  "weakly  coupled" 
cores,  it  has  been  found  that  the  transfer  function  is 
space-dependent,  e.g.,  see  reference  (Kyi  65).  Therefore, 
in  analyzing  the  kinetics  of  a weakly  coupled  system  one 
is  not  able  to  determine  a meaningful  transfer  function  on 
the  basis  of  space-independent  reactor  kinetics.  It  will 
now  be  shown  how  the  coupled  kinetics  equations  may  be 
used  to  determine  a meaningful  space-  and  energy-dependent 
transfer  function. 

Consider  the  problem  in  which  a perturbation  is 
caused  by  oscillating,  in  region  R^,  some  material  that 
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only  absorbs  neutrons  in  energy  group  E^,.  On  the  basis 
of  first  order  perturbation  theory  the  only  parameters 
directly  affected  by  the  perturbation  are  the  • 

Y ® 

The  fractional  change  in  l/e^7'n^i.s  given  by 
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(2.79) 


Equation  (2.79)  follows  directly  from  the  definition  of 


0 


given  by  Eq.  (2.57).  To  simplify  the  notation  the 


(m,n) 

k 
Y 

perturbed  parameter  will  be  written  in  the  following  way: 


l/0,fT'n)  (t)  = l/01J?'n)  (0)  {1  + e<?'n)(t)} 
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Assuming  that  small  variations  in  l/e/™'n^  produce 

^ a 

Y 

correspondingly  small  variations  in  the  neutron  population, 

we  may  expand  in  a Taylor  series  about  its  steady 

a 

state  value  and  truncate  the  series  after  the  second 
term  to  obtain: 


N.jn)  (t)  * 1 + «N^n)  (t) 

a a 


(2.82) 


Substituting  Eqs.  (2.80)  and  (2.82)  into  the  coupled 
kinetics  equations  given  by  Eqs.  (2.70)  and  (2.71), 
neglecting  terms  of  the  order  62,  Fourier  transforming 
the  resulting  expression,  and  rearranging  yields: 
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and  the  elements  of  the  column  vector  x9(w)  are  given  by 
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The  solution  of  Eq.  (2.83)  is  found  by  applying 
Cramer's  rule  (Wyl  60)  which  gives 
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where  | >{  | - 
the  vector  ^9 . 


I H I with  the  q-th  column  replaced  by 
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Summing  over  all  modes  yields  the  following 
expression  for  the  Fourier  transform  of  6N^ (t) : 
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(2.89) 


where  q = (n-1)  TK  + (a-l)I<  + l.  The  quantity  n^(w) 

a 

represents  the  output  that  would  be  measured  in  some 
region  R^  and  group  Eq  due  to  a pile  oscillator  in  region 
that  absorbs  neutrons  in  group  , . Summing  over  all 
modes  yields  the  total  input  to  region  R^,  and  group  E^,s 
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The  transfer  function  is  computed  as  the  ratio  of 
the  output  to  the  input.  In  this  case  the  transfer 
function  t7k, y t/j#y  for  an  output  in  region  R and 
group  that  is  caused  by  an  input  to  region  R^,  and 
group  E^,  is  given  byt 
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(2.91) 


where  FI  j ^ and  hO  ^ ^ are  constants  and  are  equal 
a Y ' 

to  the  corresponding  terms  in  ( } in  the  numerator  and 
denominator.  (Note:  is  pronounced  "ya"  , and  K)  is 

pronounced  "yu".) 

In  Eq.  (2.92)  the  determinant  | )-{  | that  appears  in 
the  denominator  is  the  function  of  interest.  This  function 
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describes  the  behavior  of  the  system  irrespective  of  the 
measuring  device.  In  fact,  1/ 1 )~ { | is  analogous  to  the  space- 
independent  transfer  function  derived  from  the  point-model 
kinetics  equation.  The  other  parameters  appearing  in  Eq. 
(2.91)  account  for  the  localized  nature  of  the  experimental 
arrangement  that  is  being  used  to  measure  a quantity  that 
varies  from  one  point  to  the  next.  The  behavior  of  the 
numerator  depends  on  the  location  of  the  detector  as  well  as 
the  perturbation,  and  the  behavior  of  the  denominator  depends 
on  the  location  of  the  perturbation. 


Diffusion  Theory  Form  of  the  Kinetics  Coefficients 

In  order  to  keep  the  derivation  of  the  coupled  kinetics 
equations  as  general  as  possible,  the  neutron  transport 
equation  was  used.  However,  in  a practical  problem  the 
transport  equation  is  too  difficult  to  solve  and  approximations 
must  be  made.  The  most  common  approximate  representation  of  the 
transport  problem  is  that  given  by  diffusion  theory  (Wei  58) . 
Except  for  e^m'n)  and  o<™'n)  the  mathematical  defini- 

Y r 

tions  of  the  kinetics  parameters  also  apply  to  the  diffusion 
approximation  if  we  replace  the  leakage  operator  by 
— V • D v where  D(r,E)  is  the  diffusion  coefficient.  It  is  this 
change  in  leakage  operator  that  brings  about  the  changes  in 
the  form  of  e^m'n)  and  o j™£n)  . The  diffusion  theory 
expression  for  these  parameters  is: 
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« (m,n) 
0k 
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+ 1 

( ♦ ' _ <j>  ) 

m'  v Tn;k 
Y 


{(D,7*+.V*  ) + 

m n K 


>m'  Wk  + T / dE  / 


ds, 


m(sk'E)Vsk'E)  "o  / dE  / d^k*I(®k'E)D7*n(®k'E) 


E S, 

Y k 


(2.92) 


(m,n)  _ 


<*+>  7 Vk 


i / dE  / dVm(sk'E)  VVE)  + | / dE  / dSk»^(sk,E)OT«n(Sk,E)  -n 


E S, 

Y k,  j 


E S.  . 
Y k,3 


(2.93) 


where 


Sk  - surface  of  region 

1 j ,k  = common  interface  between  regions  and  R, 
nif  = outward  drawn  normal  from  surface  S, 
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The  diffusion  coefficient  D(r,E)  is  assumed  to  be  independent 
of  r at  all  points  within  a region  but  may  vary  from  region 
to  region.  For  diffusion  theory  the  parenthesis  notation 
does  not  include  integration  over  ft. 

Conclusion 

The  theory  of  coupled  reactor  kinetics  has  been 
established  on  a firm  basis  by  deriving  the  pertinent  • 
set  of  equations  from  the  general  Boltzmann  equation. 

The  formulation  of  the  equations  from  first  principles  on 
the  basis  of  a semidirect  variational  method  allows  the 
computation  of  coupling  parameters  in  a consistent  way. 

Since  no  restrictions  have  been  placed  on  the  reactor 
configuration,  on  the  division  of  the  total  system  into 
regions,  on  the  division  of  the  spectrum  into  groups,  nor 
on  the  mathematical  representation  of  the  steady  state 
flux,  the  set  of  kinetics  equations  presented  here  is  in 
a very  suitable  form  to  study  the  dynamic  behavior  of  a 
variety  of  composite  systems. 

To  facilitate  the  utilization  of  the  coupled  formu- 
lation, the  next  chapter  is  devoted  to  the  presentation  of 
the  coupled  kinetics  parameters  in  terms  of  multigroup 
diffusion  theory. 


CHAPTER  HI 


MULTI GROUP  DIFFUSION  THEORY  FORMULATION 
OF  COUPLED  KINETICS  PARAMETERS 

Introduction 

Before  undertaking  a study  of  the  dynamic  behavior 
of  a nuclear  reactor  core,  it  is  necessary  to  study  the 
static  problem,  which  may  be  approached  analytically, 
numerically , or  some  combination  of  the  two.  There  are 
numerous  analytical  methods  of  studying  static  reactor 
physics  problems  but  these  methods  are  generally  very 
difficult  to  apply  to  composite  cores.  The  usual  approach 
is  to  approximate  the  general  transport  equation  with 
diffusion  theory  and  then  to  solve  the  diffusion  equation 
numerically.  This  method  involves  the  use  of  multigroup 
diffusion  theory.  Because  of  the  widespread  use  of  multi- 
group diffusion  theory,  a presentation  of  the  coupled 
kinetics  parameters  in  terms  of  this  formulation  should 
be  most  beneficial. 

Since  the  primary  purpose  of  this  thesis  is  to 
develop  a method  for  studying  the  kinetic  behavior  of  a 
reactor  consisting  of  coupled  cores,  a method  is  discussed 
regarding  the  numerical  determination  of  eigenmodes  for 
representing  the  flux  in  a coupled  system. 
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Multigroup  Diffusion  Equation 

The  notation  in  this  section  is  similar  to  that  used 
by  D.  C.  Bailer  in  his  description  of  the  FAIM  code  (Bal  62) . 
In  order  to  simplify  the  presentation  the  discussion  in 
this  section,  and  in  the  remainder  of  the  chapter  as  well, 
is  limited  to  a 1-dimensional  configuration  such  as  a 
slab,  sphere  or  cylinder.  The  extension  to  2-  and  3- 
dimensional  configurations  is  straightforward. 

The  multigroup  diffusion  equation  is  given  by 

i-1 

-Div2$1  (r)  + I^1(r)  = x*S (r)  + ]>  £Sfj>i  $j(r) 

j=l 

j*i  (3.1) 


for  1 < i < G.  G = number  of  groups. 


The  symbols  are  defined  as 


<t>i  (r) 


Di- 


neutron  flux  in  i-th  group  at 
position  r 

diffusion  constant  for  the  i-th 
group 

total  removal  from  the  i-th 
group  G 

D1(B2)1  + E*  + l E . . for  i <G 

a L s , l-*-] 

j*l+l 

D1 (B2 ) 1 + E 1 for  i=G 

a 

transverse  buckling  for  the  i-th 
group 

absorption  cross-section  for  the 
i-th  group 


\ 
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E 


i 

x 


r 

S (r) 


X 

(vEf)  1 


scattering  or  transfer  from 
group  i to  group  j (includes 
both  elastic  and  inelastic 
scattering) 

the  integral  of  the  fission 
spectrum  over  the  energy  range 
represented  by  group  i 

distance  measured  from  the 
origin 


fission  source  distribution 
given  by 


S (r)  = l 


(vV  i 

- — * (r) 


i=l 


X 


the  eigenvalue 

the  average  number  of  neutrons 
produced  by  a fission  in  the  i-th 
group  times  the  fission  cross- 
section  in  the  i-th  group 


Coupled  Kinetics  Parameters  in  Terms  of 
Multigroup  Diffusion  Theory- 

Before  proceeding  with  the  formulation  of  the 
coupled  kinetics  parameters  a few  remarks  should  be  made 
concerning  the  nomenclature  to  be  used.  The  terms  "broad 
region"  and  "broad  group"  will  refer  to  regions  of  the 
type  R j and  groups  of  the  type  Eq  discussed  in  Chapter  II. 
The  terms  "region"  and  "group"  refer  to  the  partitioning 
of  the  reactor  according  to  the  nuclear  properties  of 
each  section  and  to  the  groups  in  the  multigroup  structure. 
Each  broad  group  E^  consists  of  one  or  more  groups  and 
each  broad  region  R ^ may  consist  of  one  or  more  regions. 
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order  to  simplify  the  notation  the  subscripts 
and  superscripts  specifying  the  particular  eigenmode  of 
the  flux  or  adjoint  flux  are  not  shown.  The  addition  of 
this  notation  to  the  following  formulas  is  straightforward 
from  the  definitions  given  in  Chapter  II. 

E / dv  (3.: 

V UEy  Rk 

where 

$ = the  adjoint  of  x 

dV  = dr  (plane  geometry) 

= 2irr  dr  (cylindrical  geometry) 

2 

“ 4lTr  dr  (spherical  geometry) 


Let  the  weighted  production  rate  be  represented  by  G> 

JC 

where  a+y 

% 

= I I / x%+1(vEf)  V dv 

a-*y  icE  jeE  R, 

a y k (3.3) 


then  the  transfer  time  by  regeneration  is  given  by 


Ak 

a+y 


G> 


k 

a+y 


(3.4) 


The  coupled  effective  delayed  neutron  fraction  for  the 
j,-th  precursor  group  is  given  by 
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S 


l ,k 
a+y 


& 


k 

a+y 


l l 

ieE  icE 
a y 


I f^  + 1(vEf)  V dV 


(3.5) 


where 


f1 

l 


the  integral  of  the  delayed  neutron 
spectrum  for  the  i-th  precursor  group 
over  the  energy  range  represented  by 
group  i 


If  all  of  the  neutrons  emitted  by  the  l-th  precursor  group 
belong  to  some  I-th  energy  group,  then 


5il 


( *»  Kronecker  delta) 


and  Eq.  (3.5)  becomes 


$ 


l ,k 
a+y 


S>, 


k 

a-*-y 


l 

ieE 


6 • T 
ll 


l 

jeE 


4>  + I (vEf)  j<(. j dV 


(3.6) 


The  average  transfer  time  by  scattering  from  group  E 
E to  group  E is  given  by 

a y J 

<♦*'  i *>k 

Tk  = 

°^Y  III  <t>  + 1 I dV 

ieE  jeE  R.  S,j-*-i 

Y a k 


(3.7) 
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The  inverse  average  total  loss  rate  from  broad 
region  and  broad  group  E is  given  by 


U+'  k *>k 


’k 

Y 


l I + ,/iji, 

icE  R.  + » -i 

Y K R 

r=R 


where 


(3.8) 


• V 


,R 


J3- 

* 


d 
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the  outer  boundary  of  broad  region  R 

k 

the  inner  boundary  of  broad  region 

cp.,r,-{±ii£i.2i  ±>ilc)} 


2 dr 


(3.9) 


1 ; P=0 

(plane  geometry) 

2 7T  ; 

=1 

(cylindrical  geometry) 

4ir  ; 

=2 

(spherical  geometry) 

The  inverse  average  transfer  rate  into 


is  given  by 


region  R, 


(<l>+,  i.  $) 
v y k 


(k+1)  -*>k 
Y 


I (*+V( 


ieE 


r»R, 


(3.10) 
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and 


(3.11) 


In  conclusion  the  adjoint  weighted  precursor 
concentration  at  steady  state  is  given  by 


{cl'flcl)k 


(*+'  £♦>*  I 

a=l 


8 


k 

y 


l ,k 
a-»v 

Ak 

a-*-y 


(3.12) 


The  theory  as  developed  in  Chapter  II  permits  the 
description  of  the  neutron  field  in  terms  of  eigenmodes. 
In  the  next  section  this  feature  is  discussed  in  more  ‘ 


detail . 


Modal  Representation  of  the  Flux  in  a Coupled  System 

Within  the  framework  of  multigroup  diffusion  theory 
one  may  obtain  an  "exact"  solution  to  the  equations  repre- 
sented by  Eq.  (3.1)  by  solving  this  set  of  equations 
numerically.  The  accuracy  of  the  "exact"  solution  depends 
on  the  numerical  technique  that  is  used.  The  term  "exact" 
is  used  as  opposed  to  an  approximate  analytical  solution 
in  which  one  expands  the  flux  in  some  finite  series  of 
eigenfunctions  belonging  to  an  infinite  set. 
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In  the  study  of  a reactor  consisting  of  weakly 
coupled  cores  there  may  be  some  advantage  to  be  gained  by 
expanding  the  flux  representing  the  total  system  in  a 
series  of  partial  fluxes,  each  representing  some  region 
of  the  system.  One  method  of  classification  is  to  label 
the  partial  fluxes  according  to  the  region  in  which  they 
originated,  i.e., 

M 

= l ♦J(r) 

nal  (3.1 

where 

$n  2 the  flux  in  the  i-th  energy  group  due 
to  neutrons  born  in  some  n-th  region 
of  the  reactor 


The  summation  over  n includes  all  fuel  regions  of  the 
reactor.  By  definition  the  are  solutions  of  the 
following  set  of  equations: 


-DV**(r)  * £^(r) 


i-1 

y £ . , ())  3 (r) 

j=l 


xXSn(r) 


where 


(3.14) 


S (r) 
n ' ' 


G (v£.)i<^1(r) 

l — 


if  r is  the  n-th 
fuel  region 


0 otherwise 


A set  of  partial  adjoint  fluxes  is  defined  similarly; 
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The  method  presented  above  was  used  by  Avery  in  his 
"Theory  of  Coupled  Reactors"  (Ave  58a)  and  later  by 
Dougherty  in  the  development  of  his  space-time  neutron 
kinetics  equations  (Dou  61) . Also,  Stewart  has  presented 
a method  (Ste  58)  similar  to  this  one. 

The  procedure  for  determining  the  <|>*(r)  is  as  follows. 
One  first  solves  the  set  of  homogeneous  equations  expressed 
by  Eq.  (3.1)  to  determine  the  group  fluxes  4>x(r)  and  the 
eigenvalue.  He  then  uses  this  set  of  fluxes  in  the  source 
term  on  the  right  hand  side  of  Eq.  (3.14)  and  solves  the 
set  of  equations  expressed  by  Eq.  (3.14)  to  determine 
^ n ( r ) ^or  n = where  M is  the  number  of  weakly 

coupled  fuel  regions  in  the  core.  The  partial  flux  <|>*(r) 
then  represents  the  contribution  of  neutrons  originating 
in  the  n-th  fuel  region  to  the  total  flux  in  energy  group  i 
at  point  r.  To  say  that  a reactor  is  weakly  coupled  implies 
that  the  reactor  not  only  consists  of  regions  that  are  nearly 
critical  by  themselves  but  also  that  these  regions  are 
capable  of  behaving  almost  independently  during  a transient. 
It  is  for  the  latter  reason  that  it  may  be  advantageous  to 
study  weakly  coupled  reactors  in  terms  of  these  partial 
fluxes,  or  spatial  eigenmodes. 

One  very  useful  property  of  these  eigenmodes  is  that 
they  form  a complete  set  with  only  a finite  number  of  com- 
ponents in  the  sense  that  a system  {$n(r)}  is  termed 


complete  (Kap  52)  if  a function  <|>(r)  can  be  approximated 
as  closely  as  desired  by  a finite  number  of  components  4>n(r). 

Another  advantage  in  using  eigenmodes  that  are 
representative  of  a particular  region  of  the  reactor  is 
that  the  number  of  coupled  kinetics  equations  for  a given 
problem  may  be  greatly  reduced.  This  results  from  the 

restriction  that  the  precursors  do  not  diffuse  which  requires 
that 


(CA  #m' f tCi ,n)k' 


6 (cT  c \ 

m,nv  £,n'rlc*,n;k 


5 U+  i $ ) Y 
m,n ' m v n'k  1 


y ct=l  A 


, (m,n) 
'l,k 

a-»Y 

(m,n) 

k 

a-*-y 


(3.15) 


The  quantity  c£,n  gives  the  number  of  neutrons  emitted 
into  the  n-th  spatial  mode  by  the  £-th  precursor  group. 

The  quantity  c^n  i3  the  adjoint  to  c£  ^ and  both  quantities 
are  non-zero  only  within  the  fuel  region  specified  by  the 
subscript  n.  Furthermore,  S^n)  is  non-zero  only  if 

the  n-th  fuel  region  is  within  broad  region  R^.  in  the 
general  case  discussed  in  Chapter  II  it  was  necessary  to 
solve  up  to  MrK ( 1+D)  coupled  first  order  differential 
equations.  In  this  case  it  is  only  necessary  to  solve  up 
to  TK(M+D)  coupled  first  order  differential  equations. 
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One  disadvantage  of  the  method  discussed  above  is 
that  most  multigroup  diffusion  theory  codes  are  not 
programmed  to  perform  the  calculations  necessary  for 
determining  the  partial  fluxes  and  their  ad joints.  However, 
this  problem  is  not  a real  hindrance  since  with  a few 
modifications  some  existing  code  can  be  made  up  to  perform 
the  additional  operations.  Such  modifications  have  been 
made  in  the  FAIM  code  (Bal  62)  and  the  modified  code  has 
been  given  the  name  SAIM  (Coc  65a)  . 

The  SAIM  code  computes  values  for  the  partial  fluxes 
and  partial  adjoints  and  punches  them  on  cards  along  with 
the  corresponding  cross-sections.  This  punched  card  output 
is  then  read  directly  into  a numerical  integration  code 
called  C0UPKIN  which  computes  the  coupled  kinetics  parameters 
according  to  Eq.  (3.2)  - (3.11). 

The  C0UPKIN  code  is  dimensioned  to  compute  the 
coupling  parameters  for  a system  having  up  to  twenty  broad 
regions  with  the  energy-dependence  expressed  by  up  to  eighteen 
broad  groups.  A complete  description  of  the  code  is  given 
in  reference  (Coc  65b)  . 

Conclusions 

Formulas  have  been  given  for  computing  the  coupled 
kinetics  parameters  in  terms  of  multigroup  diffusion  theory, 
and  a method  has  been  suggested  for  studying  weakly  coupled 
systems  by  expanding  the  flux  in  spatial  eigenmodes. 
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In  the  next  chapter  the  coupled  kinetics  formulas  are 
applied  to  a kinetics  study  of  a simplified  model  of  an 
Argonaut-type  reactor  with  a two-slab  loading. 


CHAPTER  IV 


SPACE-TIME  ANALYSIS  OF  THE  PULSED  NEUTRON-BEAM 
EXPERIMENT  IN  A REACTOR  WITH  WEAKLY 
COUPLED  FUEL  REGIONS 

Introduction 

One  type  of  experiment  that  is  well  suited  to  the 

\ 

purpose  of  determining  certain  fundamental  reactor  parameters 
is  based  on  the  pulsed  neutron-beam  technique  (SjO  56) . 

This  experiment  consists  essentially  of  emitting  a very 
short  burst  of  neutrons  into  the  reactor  and  measuring  the 
decay  constant  of  the  most  slowly  decaying  mode.  Then  on 
the  basis  of  a 1-point  reactor  kinetics  model  this  measured 
decay  constant  is  related  to  certain  of  the  reactor  parameters. 
Typically,  one  uses  this  technique  to  determine  the  ratio  of 
the  effective  delayed  neutron  fraction,  B,  to  the  mean  gener- 
ation time,  A,  at  delayed  critical,  and  then  by  assuming  that 
S/A  is  independent  of  the  reactivity,  it  is  possible  to  make 
reactivity  measurements.  In  the  study  of  reactors  with  a 
single  fuel  region  analysis  of  the  experimental  data  on 
the  basis  of  1- point-model  kinetic  theory  has  been  shown  to 
give  excellent  agreement  with  theoretical  computations, 
for  example  see  reference  (Per  60) . However,  when  theoretical 
computations  were  compared  with  experimental  measurements 
made  on  an  Argonaut-type  reactor  with  a 2-slab  loading, 
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the  comparison  was  so  poor  that  questions  arose  regarding 
both  the  experimental  measurements  and  the  theoretical 
computations.  Among  these  questions  were  the  following; 

(a)  Will  1-point-model  reactor  kinetics  correctly  predict 
the  prompt  neutron  decay  constant  for  a reactor  with  weakly 
coupled  fuel  regions?  (b)  Will  the  coupled  kinetics  for- 
mulation give  a better  prediction?  (c)  How  well  do  these 
point-model  approximations  agree  with  the  true  value?  (d) 

What  effect  does  the  graphite  in  the  region  between  the  fuel 
slabs  have  on  the  value  of  the  decay  constant?  (e)  Is  the 
experimenter  actually  measuring  the  prompt  neutron  decay 
constant , or  is  he  measuring  the  die-away  of  the  pulsed 
neutron  beam  in  the  graphite? 

An  answer  to  these  questions  was  sought  through  the 
study  of  a very  simplified  model  of  a reactor  with  two 
weakly  coupled  fuel  regions.  This  model  is  shown  schemat- 
ically in  Figure  1.  One  important  feature  that  has  been 
retained  in  this  model  is  the  separation  of  the  two  fuel 
regions,  1 and  3 , by  a reflector  region,  2.  Three  methods 
of  analysis  will  be  presented.  The  first  method  will  be 
considered  an  exact  representation  of  this  model  which  is 
based  on  the  time-  and  space-dependent  monoenergetic  diffusion 
theory  model.  The  second  method  will  utilize  the  coupled 
kinetics  formulation  and  the  third  method  will  present  the 
1-point-model  representation.  Numerical  computations  of 
the  prompt  neutron  decay  constant  as  predicted  by  the 
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Vacuum 


Figure  1. 
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Reactor  With  Two  Coupled  Fuel  Regions 
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three  methods  will  be  compared  in  order  to  answer  questions 
(a),  (b) , and  (c) . Also,  in  the  process  of  obtaining  these 
numerical  values  the  answers  to  questions  (d)  and  (e)  will 
be  revealed. 

Space-  and  Time-Dependent  Reactor  Kinetics 
Representation  " , 

Consider  the  following  problem.  A reactor  is  delayed 
critical  and  operating  at  steady  state.  At  time  t=0  a high 
intensity  neutron  beam  of  short  duration  is  introduced.  The 
neutrons  from  this  external  source  cause  fissions  which  pro- 
duce more  neutrons,  both  prompt  and  delayed.  The  net  result, 
after  the  delayed  neutrons  begin  to  appear,  is  an  increase 
in  the  steady-state  flux.  However,  during  the  short  time 
interval  between  the  introduction  of  the  neutron  beam  and 
the  appearance  of  the  additional  delayed  neutrons,  there  is 
a measurable  decrease  in  the  neutron  population  as  a function 
of  time  which  exponentially  approaches  a steady-state  dis- 
tribution. At  first,  this  decrease  in  the  neutron  popula- 
tion is  due  primarily  to  the  loss  of  the  injected  neutrons 
by  absorption  and  leakage,  but  it  will  be  shown  that  the 
duration  of  this  effect  is  brief  compared  to  the  duration 
of  the  phenomena  in  which  we  are  interested.  The  phenomenon 
of  interest  is  the  rate  at  which  the  prompt  neutron  popula- 
tion, or  flux,  that  resulted  from  fissions  caused  by  the 
injected  neutrons,  dies  away.  The  die-away  occurs  because 
the  prompt  neutrons  are  unable  to  sustain  a chain  reaction 
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in  a delayed  critical  reactor.  These  prompt  neutrons  will 
be  referred  to  as  "excess  prompt  neutrons"  in  the  discussion 
that  follows: 

Following  the  injection  of  the  neutron  beam  at  time 
t=0,  the  flux  in  the  reactor  consists  of  a steady-state 
portion  $(x)  which  also  existed  at  t<0  and  a time  varying 
portion  made  up  of  the  excess  prompt  neutron  flux  i|i(x,t) 
and  the  source  flux  s(x,t|0,0).  The  function  s(x,t|0,0) 
gives  the  flux  of  source  neutrons  at  x and  t due  to  injec- 
tion of  a neutron  beam  at  x=0  and  t»0. 

For  the  delayed  critical  reactor  shown  in  Figure  1 
the  steady-state  flux  must  satisfy  the  following  set  of 
equations  s 


V 2 4>  ( x ) 

+ B2<Kx)  = 0 

b< | x| <a 

(4.1) 

V2*(x) 

- k 2 <j>  (x)  = 0 

0< | x | <b 

(4.2) 

where 

B2  = 
m 

U£f  ' £ RF 

(4.3) 

°F 

2 

r a 

erg 

dg 

(4.4) 

The  subscripts 

F and  G refer  to 

the  fuel  regions  and 

the 

graphite  reflector  between  the  fuel  regions,  respectively. 
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The  removal  cross-section  includes  losses  due  to  absorption 
and  transverse  leakage. 

The  following  boundary  conditions  must  be  satisfied 
at  all  times: 


♦ ( ta)  - 0 

(4.5) 

e1^ra0  * (+{b+e  })  , cllmQ  4>  (±  { b-e  } ) 

(4.6) 

e^m0  DpV$  (x)  = eHm0  DGV<fr(x) 

x=  ±(b+e)  x= 

± (b-e ) 

(4.7) 

The  solution  of  Eqs.  (4.1)  and  (4.2)  is 
forward  and  is  given  by: 

straight- 

♦(*)  " Al3in  Bm(a-|x|)  b<|x|<a 

(4.8) 

(x)  = A2cosh  k|x|  0<  | x | <b 

(4.9) 

where  the  coefficients  A^  and  A2  satisfy  the  following  set 
of  homogeneous  equations: 

qX  = 0 

(4.10) 

A = [ A^  A2  ] 

(4.11) 

69 


sin  B (a-b) 
m 


-cosh  k b 


Q 


DPB  cos  B (a-b) 

r m m 


DgK  sinh  < b 


(4.12) 


A non-trivial  solution  exists  for  Eq.  (4.9)  only  if  the 
determinant  of  Q is  identically  zero.  This  is  accomplished 
by  adjusting  the  macroscopic  cross-sections.  Let  us  proceed 
to  the  solution  of  excess  prompt  neutron  flux. 

The  flux  represented  by  <>>(x,t)  is  given  by  the 
solution  to  the  following  set  of  equations: 


1 3 

— <J»(x,t) 

V 3t 


(1-6) vE 


f 


S(x,t|0,0) 


D 


F 


b< | x | <a 


(4.13) 


V2\|>  (x,t)  -<2iHx,t)  = -i—  — it/  (x,t) 

°GV  3t 


0 < | x | <b 


(4.14) 


where 


(1-8) vE 


f 


RF 


B 


2 


D 


F 


(4.15) 
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Laplace  transforming  Eqs.  (4.13)  and  (4.14)  yields: 

2—  -2—  (l“6)vE» 

V t(x,s)  + B iMx,s)  “ - S (x,s|  0,0) 

D 

F 

b< | x | <a  (4.16) 

72^(x,s)  -k2Hx#s)  = 0 0<  | x | <b  (4.17) 


where 


-2  2 
B » B - 


(l-g)vEf  - Erf  - s/v 


°pv 


D, 


(4.18) 


< 2 a * 2 x 5 ERG  + S/V 

"T*  ■ ' S5  ■ ii  . — i ■ ii 

V °G 


(4.19) 


00  -st 

?(x,s)  = / ip(x,t)  e dt 

o 


-st 


S(x,s|0,0)  = / S(x,t|0,0)  e dt 


Equations  (4.16)  and  (4.17)  must  satisfy  the  boundary 
conditions  given  by  Eqs.  (4.5)-  (4.7).  Application  of  these 
boundary  conditions  yields  the  following  solution  for  the 
Laplace  transform  of  the  excess  prompt  neutron  flux: 
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*(x,s)  * C-^s)  sin  B ( a-  | x | ) + * (xfs) 


(4.20) 


\i>(x,s)  = C2(s)  cosh  k\x\ 


(4.21) 


where  ?p(x,s)  is  the  particular  solution  that  results 
from  the  source  term  on  the  right  hand  side  of  Eq.  (4.16). 
The  solution  5»p(x,s)  represents  the  flux  due  to  those 
neutrons  that  are  produced  directly  by  fissions  caused 
by  the  injected  source  neutrons.  The  coefficients  of  the 
complementary  solution  are  given  by  the  solution  to  the 
following  set  of  equations: 

Z (s)  C (s)  =*  W(s)  (4  2 


where 


£(s)  - [C1 ( s)  C2 (g) 3 


(4.23) 


sin  B(a-b) 


-cosh  k b 


Z (s) 


DpB  cos  B(a-b)  Dg<  sinh  < b 


(4;24) 


P 


W(s) 


DF7fp(b+,s) 


(4.25) 
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Using  Cramer's  rule  (Kap  52)  to  solve  Eq.  (4.25)  yields: 


C.(s) 


IVs) 1 

| Z (s)  | 


i-1,2 


(4.26) 


where  |z^(s) | is  the  determinant  of  the  matrix  Z(s)  with 
the  i-th  column  replaced  by  W(s) . 

The  time-  and  space-dependent  solution  of  Eqs.  (4.13) 
and  (4.14)  is  given  by  the  inverse  Laplace  transform  of 
£(x, s)  : 


1 

* (x,t)  « 

2ni 


o+iB  _ st 

/ ?(x,s)  e ds 


(4.27) 


Substituting  Eqs.  (4.20)  and  (4.21)  into  Eq.  (4.27)  and 
using  the  expression  for  C^s)  and  C2(s)  given  by  Eq.  (4.26) 
yields : 


1 

*(x,t)  - 

2iri 


o+i 


o-i 


00 


| Z ! ( s ) I 
| z (s)  I 


sin  B ( a—  J x| ) 


st 

e ds 


+ 


1 


2-ni 


o+i 
/ " 
o-i 


st 

ifp  (x , s ) e ds 


b < | x | <a 


(4.28) 
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♦ (x,t) 


St 

cosh  ic  I x I e ds 


CO 


(4.29) 


The  integral  over  *p(x,s)  gives  the  contribution 
from  fissions  caused  directly  by  the  source  neutrons  which 
die  away  rapidly  compared  to  the  asymptotic  mode  of  the 
excess  prompt  neutron  flux  which  is  given  by  the  residue 
at  the  smallest  zero  of  |z(s) | . It  is  the  decay  constant 
of  the  asymptotic  mode  that  is  measured  in  the  pulsed 
neutron-beam  experiment.  The  asymptotic  solutions  of  Eqs. 
(4.28)  and  (4.29)  are  given  by: 


asym 


(x,t)  « 


sin  Bq ( a- | x | ) e 


d 


— | Z (s)  | 
ds  s=s 


b<  | x | <a 


(4.30) 


(X,t)  = 


asym 


cosh  j(0 


I x | e 


ds 


s=s 


0 


0< | x | <b 


(4.31) 


5 


0 


where 
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K 

O 


The  value  of  sq  is  given  by  the  smallest  root  of  the 
following  transcendental  equation: 


tan  B(a-b) 


dfb 


coth  < b 


°G  ^ 


(4.32) 


Equation  (4.32)  was  obtained  by  setting  |Z(s)|  =0. 

In  order  to  continue  this  study  let  us  assume 
some  realistic  values  for  the  system  parameters.  The 
University  of  Florida  Training  Reactor  (Dun  58)  will  be 
used  as  a guide.  In  this  case  the  following  values  are 
assumed : 


a 

83 

30  cm 

erf 

as 

0.10  cm  ^ 

b 

ss 

15  cm 

df 

ss 

0.25  cm 

K 

83 

0.04  cm  ^ 

8 

ss 

0.0064 

dg 

SS 

1.00  cm 

V 

S3 

2.2  x 10^cm/sec 

(4.33) 


The  reflector  region  between  the  fuel  slabs  in  the  UFTR 
is  graphite  which  has  a value  of  0.02  cm-1  for  the  quantity 
represented  by  < • However,  to  make  the  values  used  in 
this  model  more  realistic  the  value  of  < has  been  increased 
to  account  for  transverse  leakage.  Likewise,  the  value  of 
ERF  is  greater  than  the  value  for  the  thermal  absorption 
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cross-section  in  the  fuel  region  of  the  UFTR  in  order  to 

account  for  absorption  while  slowing  down  and  for  trans- 
verse leakage. 

When  these  values  are  inserted  into  the  determinant 

of  Q given  by  Eq.  (4.12),  the  smallest  value  of  B which 

m 

satisfies  the  critical  condition  that  |q|  = 0 is 

Bm  = °*1412  cm~1  (4.34) 

which  requires  that 

= 0.10498  cm-1  (4.35) 

and 

— - 0.7205 

A1 

Inserting  the  values  given  by  Eqs.  (4.33)  and  (4.35)  into 
Eqs.  (4,18)  and  (4.19)  yields  the  following  expressions 

for  B and  < : 

% 

B = 0.00426  /~94  8-s~ 

< - 0.00213  /352+s  (4.36) 

The  root  of  |z(s)|  with  the  smallest  magnitude  is  at 

so  = ”73*5  sec-1  (4.37) 

Thus,  the  decay  constant  of  the  most  slowly  decaying  mode 
of  the  excess  prompt  neutron  flux  is  equal  to  -73.5  sec”1. 
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In  the  next  section  the  coupled  kinetics  formulation 
presented  in  Chapter  II  is  used  to  analyze  this  problem 
in  terms  of  a 3-point  model. 

Coupled  Point-Model  Representation 

In  terms  of  the  coupled  point-model  representation 
the  equations  describing  the  kinetics  of  the  flux  in  the 
reactor  shown  in  Figure  1 are  as  follows: 


d 

~ N.  (t) 
dt  1 


1-6 

[ + 


1 


1 1 

]N1(t)  -I-  - N2(t) 

9p  a 

G+F  (4.38) 


d 

— N2(t)  » 

dt 


(t) 


1 1 
+ - Nx(t)  + - N 3 ( t ) 

°F-*-G  •>  aF-*-G 

(4.39) 


d 

— N,(t) 
dt  J 


1 1 

— + — ] N (t) 
6„  t _ J 


1 

+ — 

°G-+-F 


n2  (t) 


(4.40) 


whe  re 


the  generation  time  of  the  neutrons 
in  region  1 or  3 


the  loss  rate  of  neutrons  due  to 
leakage  and  absorption.  Subscript 
F refers  to  the  fuel  regions  and 
subscript  G refers  to  the  reflector 
region 
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1 


°G->F 

1 


°OF 


the  transfer  rate  of  neutrons  from 
the  reflector  to  the  fuel 


the  transfer  rate  of  neutrons  from 
the  fuel  to  the  reflector 


For  simplicity  only  one  space  mode  is  considered. 

The  formulas  for  evaluating  the  parameters  follow 
directly  from  the  general  theory  presented  in  Chapter  II. 
For  clarity  they  are  presented  here  in  terms  of  this 
simple  example. 


4f,  I *>F 


A sl 

i / sin1 2B  (a-x)dx 
T ' m 


v 


(4.41) 


(*+,  i *>G 


2A 


2 b 


2 

/ cosli  < x dx 


(4.42) 


1 U ,vEf<fr)p 

af  U+,  \ *)F 


= vvE 


(4.43) 


♦ (b) j_(b) 


F-^G 


(/,  i *)G 


(4.44) 


1 *T(b)j+(b) 


G->F 


(♦f,  I *)F 


(4.45) 


78 


1 1 

[ — ] = VE 

6F  TF 


(V4,+  #Dp74»)F  *+(b)jJb) 


RF 


<*+'  7 *>F 


* ' V ^ ^ F 


(4.46) 


1 1 

[ 3 = vE 

6G  TG 


RM 


(V^+,Dg^)g  2 <|> + (b)  j+  (b) 

+ . 

(♦t  I f).  ($+  I $) 

»vy/G  ' v G 


(4.47) 


4>  (x)  = <j>  (x) 


(4.48) 


A^cosh  < b A»D_  < 3inh  < b 

j_(b)  * — t JL G 

+ 4 2 

A.sinB  (a-b)  A,D_B  cosB  (a-b) 

„ m ~ i t m m 

4 2 

(4.49) 

Expressions  for  <j>  (x)  are  given  by  Eqs.  (4.11)  and 
(4.12) . The  parenthesis  notation  indicates  integration  over 
the  region  specified  by  the  subscript. 

Inserting  the  values  given  by  Eqs.  (4. 36) - (4 . 38) 
into  Eqs.  ( 4 . 41 ) - ( 4 . 49)  yields  the  following  values  for 
the  kinetics  parameters: 
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<*+.  £ *>p  - 4.125  X 10*5  A*  (4.50) 

% 

U+,  ^ *>G  = 1,540  X 10~4  A2 


1-8 

af 

= 2.2948  x 104  sec"1 

1 

°F-G 

- 2.378  x 103  sec"1 

1 

°G-F 

= 4.224  x 103  sec"1 

1 

c 

°F 

■ — ] = 2.770  x 104  sec"1 

tf 

1 

c 

9g 

■ — ] = 4.752  x 104  sec"1 

tg 

Laplace  transforming  Eqs.  (4 . 38) - (4 . 40)  yields 

H(s)  n (s)  = | (4.51) 


-► 

n *= 

r , (4.52) 

L (s)  r>2  (s)  (s)  ] 

t = 

Cl  1 1]  (4.53) 
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(s)  = 


, d-8)  1 1 

{s  + __} 


tf  6f 


F-+G 


G+F 


1 1 

{s + — } 


TG  eG 


F+G 


G-+F 


{s  - 


(1-8) 


1 1 
— + — } 
TF  0F 


n.(s)  = / N (t)  e 

o 


st 


dt 


j = 1,2,3 


(4,54) 


(4.55) 


Application  of  Cramer's  rule  yields  an  expression 


for  the  n j (s) : 


IK  j<s>  I 

1 j ( s)  = 2 

I H (s)  | 


(4.56) 


where  |Hj (s) | is  the  determinant  of  H(s)  with  the  j-th 
column  replaced  by  t.  The  solution  of  Eqs.  ( 4 . 38 ) - ( 4 . 40) 
is  given  by  the  inverse  Laplace  transform  of  the  nj(s)j 


N.  (t) 


1 

2 u i 


o+i“ 

/ 

o-i“ 


I K -i  < s>  I st 

e ds 

I H(s)  | 


(4.57) 


81 


The  solutions  of  the  (t)  are  given  by  the  sum  of 
three  exponential  terms  determined  by  the  residue  at  the 
three  zeros  of  |H(s)|.  The  smallest  root  of  |H  (s) | is 

s ■ 74.8  sec"*'*' 

o 

Thus , the  value  for  the  prompt  neutron  decay  constant  pre- 
dicted by  the  coupled  kinetics  formulation  is  1.77  per  cent 
higher  than  the  value  predicted  by  the  exact  method. 

Let  us  now  proceed  to  an  analysis  of  the  problem 
in  terms  of  the  1-point-model  kinetics  equation. 

1-Point-Model  Representation 

The  1-point-model  kinetics  equation  for  a delayed 
critical  reactor  is  given  by 

d 3 

— N ( t)  N(t) 

dt  A (4.58) 

where  X is  the  mean  generation  time  averaged  over  the 
entire  reactor.  Note  that  delayed  neutrons  have  been 
neglected  as  in  the  first  two  analyses. 

The  solution  of  Eq.  (4.63)  is 

-Bt 

N(t)  = N(0)  e A (4.59) 

The  formula  for  computing  - i is  as  follows: 
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BU  'vEf4’)  REACTOR 


( * + . — <k ) 

* ' v v ' REACTOR 


26(^+,vZf^)F 

at  — . 

2(*  + , i Op  + U+,  i *)G 

$ 

(*f,  i *)G 

1 + 

2U+,  i *)F 

147.8 

S — - - - 

0.492 
1 + 

0.508 

= - 75.1  sec”1  (4.60) 

We  see  from  Eq.  (4.60)  that  the  region  between  the  fuel 
slabs  has  the  effect  of  decreasing  the  magnitude  of  the 
prompt  neutron  decay  constant  and  increasing  the  mean 
generation  time  of  the  reactor.  Note  that  approximately 
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half  the  weighted  neutron  population  is  in  region  2. 

This  doubles  the  mean  generation  time  and  as  a result 
slows  down  the  system  response. 

Comparing  the  value  of  the  prompt  neutron  decay 
constant  predicted  by  the  1-point-model  kinetics  equation 
with  the  exact  value  indicates  an  error  of  2.2  per  cent. 

Discussion  of  the  Results 

The  prompt  neutron  decay  constant  as  predicted  by 
the  three  models  considered  in  this  study  is  summarized  in 
the  following  table: 


TABLE  1 


SUMMARY  OF  RESULTS 


Model 

Decay  Constant 

Error 

Exact  , 

-73.5  sec”*- 

- 

3-point 

-74.8  3ec"^ 

1.77% 

1-point 

-75.1  sec"'*' 

2.20% 

The  agreement  between  the  three  methods  of  computation  is 
remarkable  and  speaks  well  for  the  variational  principle 
upon  which  the  point-model  representations  are  based.  > 
Recall  from  Chapter  II  that  both  the  1-point-model  and 
the  3— point— model  are  based  on  the  general  theory  of 
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coupled  point-model  reactor  kinetics  developed  in  Chapter  II. 
This  study  indicates  that,  at  least  for  this  1-group  approxi- 
mation, spatial  effects  are  not  an  important  consideration  in 
predicting  the  asymptotic  behavior  of  the  prompt  neutrons. 
Since  the  1-point-model  does  very  accurately  predict  the 
prompt  neutron  decay  constant,  then  Rossi-a  measurements  are 
meaningful  even  for  this  coupled  system  and  one  may  deduce 
from  measurements  of  the  Rossi-a  the  quantities  p and  B/A 
for  the  system.  It  must  be  emphasized,  however,  that  broad 
generalizations  should  not  be  made  on  so  simple  a model. 

This  model  is  indeed  an  over-simplification  of  the  true 
system.  However,  in  the  next  chapter  the  results  of  an 
18-group  computation  performed  on  the  University  of  Florida 
Training  Reactor  are  presented,  and  it  is  shown  there,  also, 
that  an  "exact"  numerical  computation  and  the  1-point  model 
are  in  close  agreement. 

This  study  has  indicated  the  importance  of  the  graphite 
reflector  in  increasing  the  mean  generation  time.  As  shown 
in  Eq.  (4.60)  approximately  one-half  of  the  weighted  neutron 
population  was  in  the  reflector  between  the  fuel  regions. 

This  doubled  the  mean  generation  time.  This  effect  requires 
that  a certain  precaution  be  taken  when  an  equivalent  1- 
dimensional  model  is  used  to  represent  a core  of  this  type. 

In  a 1-dimensional  representation  the  transverse  leakage 
is  assumed  to  be  proportional  to  the  flux.  The  constant 
of  proportionality  is  taken  to  be  the  product  of  the 
diffusion  coefficient  times  the  transverse  geometric 


85 


buckling  (Gla  52) . Now,  the  transverse  buckling  is  not 
easily  determined  for  a reactor  of  this  type  and  this 
creates  a problem  when  one  attempts  a 1-dimensional  compu- 
tation. If  a single  value  is  assumed  for  the  transverse 
buckling  for  all  groups  in  a multigroup  calculation,  this 
single. value  will  over-estimate  the  removal  cross-section 
for  the  thermal  group  because  the  transverse  buckling  for  • 
the  thermal  group  is  smaller  than  this  "average"  value  that 
is  assumed  for  all  groups.  Recall  from  Chapter  III  that 

the  total  removal  cross-section  iT  for  the  thermal  group 
is  given  by 

ET  “ Ea  + Db2 

2 

In  the  fuel  the  transverse  leakage  term  DB  is  small 
compared  to  the  large  thermal  absorption  cross-section  l . 

cl 

As  a result,  the  total  removal  cross-section  in  the  fuel 
is  not  strongly  affected  by  an  error  in  the  transverse 
buckling.  However,  in  the  graphite  the  thermal  absorption 
cross-section  is  very  small  and  the  value  of  the  total 
removal  cross-section  is  strongly  dependent  on  the  trans- 
verse leakage  term.  The  significance  of  this  source  of 
error  is  demonstrated  in  the  next  chapter. 

One  additional  problem  that  must  be  considered  is 

whether  or  not  the  asymptotic  decay  constant  is  really  due 

% 

to  the  prompt  neutron  die-away.  It  is  conceivable  that  one 
might  actually  be  measuring  the  die-away  of  the  pulsed 
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neutron-beam  because  of  the  low  thermal  absorption 
cross-section  of  graphite.  The  values  in  Table  2 taken 
from  reference  (Per  65)  will  be  used  to  investigate  this 
. possibility. 


TABLE  2 

THERMAL  NEUTRON  PARAMETERS  AS  DETERMINED  FROM 
VARIOUS  DIE-AWAY  EXPERIMENTS  IN  GRAPHITE 


Investigator 


aQ  ( sec"'*') 


D x 10 
o 

, 2 -1. 

(cm  sec  ) 


Lalande  (Lai 

61) 

79.5 

t 

3 

2.07 

± 

.04 

Starr  & Price 

(Sta 

62) 

Type  GBF 

74.6 

t 

.06 

2.14 

± 

.01 

Type  AA 

68.1 

t 

.08 

2.19 

± 

.03 

Klose,  et  al. 

(Klo 

62) 

88.3 

i 

1.2 

2.13 

t 

.02 

Sagot,  et  al. 

(Sag 

64) 

74.1 

i 

3 

2.19 

± 

.03 

The  decay  constant  of  the  source  is  given  by 

“SOURCE  = " “o  ' DoB  + 0(B  * (4.61) 

2 

where  aQ  and  Dq  are  given  in  Table  2 and  B is  the  geometric 
buckling.  Consider  the  graphite  region  between  the  fuel 
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slabs.  The  geometric  buckling  of  this  region  can  be  very 
well  approximated ( by 


2b 


where  2b  is  the  distance  between  the  fuel  slabs.  The  value 
of  2b  for  the  model  used  in  this  chapter  is  30  cm  (for  the 


UFTR  2b  = 30.48  cm)  which  gives  a value  of 
2 ir  2 

B » ( — ) = 0.011  cnt 

30 


(4.62) 


From  Table  2 it  is  seen  that 
ao  * ^ sec  ^ 


(4.63) 


Substituting  the  values  given  by  Eqs.  (4 . 62) - (4 . 64)  into 
Eq.  (4.61)  yields 


“source  “ 


- 2385  sec 


Therefore , the  pulsed  neutron-beam  dies  away  quickly  and 
should  not  contaminate  the  asymptotic  die-away  of  the 
prompt  neutrons  even  if  the  measurements  are  made  with  a 
detector  located  in  the  graphite  at  the  center  of  the 


reactor. 
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Conclusions 

A 1-group,  1-dimensional  model  of  a reactor  with 
two  weakly  coupled  fuel  regions  has  been  studied.  The  results 
show  that  the  1-point-model  kinetics  equation  gives  a 
reasonably  accurate  estimate  of  the  asymptotic  decay 
constant  in  a prompt  neutron  die-away  experiment.  This  result 
tends  to  indicate  that  spatial  effects  are  not  so  great  as 
to  prohibit  the  use  of  a 1-point-model  reactor  kinetics'  for 
analyzing  the  data' from  a pulsed  neutron-beam  experiment. 
Therefore,  one  may  deduce  meaningful  system-averaged  values 
for  b/A  and  the  reactivity  from  measurements  of  the  asymptotic 
decay  constant. 

The  fact  that  the  decay  constant  predicted  by  the 
3-point  model  was  closer  to  the  exact  value  than  was  the 
value  predicted  by  the  1-point  model  indicates  that  the 
coupled  formulation  converges  toward  the  exact  representation 
as  the  number  of  nodes  is  increased.  This  result  gives  con- 
fidence in  the  formulas  used  to  compute  the  coupling  parameters. 
Actually,  using  the  coupled  kinetics  equations  to  determine  the 
asymptotic  decay  constant  introduces  unnecessary  complexity 
to  the  computation.  The  coupled  kinetics  equations  are 
better  suited  for  studying  space-  and  energy-dependent  effects. 
In  fact  the  virtue  of  using  this  method  is  that  it  provides  a 
means  of  extending  static  space-  and  energy-dependent  flux 
computations  into  the  time  domain  by  solving  a set  of  first 
order  differential  equations. 


CHAPTER  V 


A DISCUSSION  OF  COMPUTER  CALCULATIONS 
AND  A SUMMARY  OF  THE  RESULTS 


Introduction 

In  this  chapter  the  results  of  a numerical  study  of 
the  neutron  flux  in  the  University  of  Florida  Training 
Reactor  are  presented.  Both  a 4-group , 2-dimensional  and 
an  18-group , 1-dimensional  computation  of  the  flux  have 
been  performed.  The  results  of  the  computations  are  used 
to  determine  the  effective  delayed  neutron  fraction,  the 
mean  generation  time  and  the  kinetics  parameters  for 
several  coupled  point  models.  A method  is  also  presented 
for  computing  the  asymptotic  decay  constant  that  is 
measured  in  a prompt  neutron  die-away  experiment.  This 
method  is  considered  "exact"  in  the  sense  that  it  does  not 
require  the  assumption  of  separability  that  is  required  in 
point-model  kinetics  calculations.  The  decay  constant  pre- 
dicted by  the  exact  method  is  compared  with  the  decay  con- 
stant predicted  by  the  1-point  model,  as  well  as  by 
several  coupled  multiple-point  models.  This  comparison  is 
made  in  order  to  check  the  accuracy  of  the  theory  presented 
in  Chapter  II.  To  check  the  accuracy  of  the  18-group,  1- 
dimensional  numerical  computations,  the  calculated  value  of 
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the  decay  constant  is  compared  with  the  values  measured  by 
a pulsed  neutron-beam  experiment  and  by  a random  noise 
experiment. 

In  addition  to  numerical  results  a brief  discussion 
of  the  computer  codes  that  were  used  is  also  presented.  In 
the  concluding  section  a procedure  is  suggested  for  deter- 
mining a region-dependent  transfer  function. 

Exact  Determination  of  the  Prompt  Decay  Constant 

Since  it  is  the  smallest  decay  constant  that  is  of 
interest  in  a prompt  neutron  die-away  experiment,  it  is 
possible  to  use  a variational  method  to  determine  this 
parameter.  The  term  "exact"  is  used  in  the  sense  that 
separability  of  time  and  space  is  not  assumed.  As  far  as 
the  "accuracy"  of  the  method  is  concerned,  it  is  only  as 
good  as  the  model  used  to  represent  the  UFTR,  namely  multi- 
group diffusion  theory. 

Following  the  injection  of  a burst  of  neutrons  into 
the  reactor,  the  excess  prompt  neutron  flux  is  described  by: 

1 3 + 

- — i/i  (r  ,E,t)  = H*(r,E,t) 

v 3t  (5.1) 

The  "excess"  prompt  neutrons  are  those  neutrons  that  resulted 
from  fissions  caused  by  the  pulsed  neutron-beam  source.  The 
operator  H is  given  by 

H*(r,E,t)  = d-8)  xp  (E)  / dE'vZf(r,E,)<l»(r,E'  ,t)  + 

E' 
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V«D(r,E)v*(r,E,t)  - ET(r,E)*(r#E,t)  + 

/ dE'E (?,E'-*E)  * (r,E't) 

E'  (5.2) 

where 

D 

(l-B)Xp(E)  = X(E)  - l fi(E)Bi 

i=l 

Xp(E)  = prompt  neutron  fission  spectrum 

X (E)  = fission  spectrum  for  all  neutrons 

resulting  from  a fission 

The  other  symbols  have  the  usual  meanings  and  have  been 
defined  in  preceding  chapters. 

Laplace  transforming  Eq.  (5.1)  yields: 

i 

HtjJ  (r,E,s)  = i ij/(r,E,s)  (5  3) 

In  order  to  minimize  the  error  in  determining  s from 
Eq.  (5.3)  a variational  method  will  be  used.  Multiplying 
Eq.  (5.3)  by  a weighting  function  £+(r,E,s)  and  integrating 
over  all  values  of  E and  r yields: 

/ dE  / dr  ?+  (r,E,s)H$(r,E,s)  = s / dE  / dr  JJJ  +( r ,E , s)  i-  £(r,E,s) 


(5.4) 
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Taking  variations  with  respect  to  s,  and  £ and  requiring 
6s  to  vanish  yields: 


/ dE  / dr  [6*+{h£  - £}  + 6iI{H+^+  - ^ } ] = 0 


(5.5) 


Tt  • 


where  the  operator  H is  the  adjoint  of  H and  is  defined 
by  (Sel  58) : 


/ gH+f  dx 

all  x 


/ fHg  dx 

all  x 

.rt 


Applying  the  usual  argument  that  5 IJJ  and  6\ji  in  general  are 
non-zero  everywhere  except  at  the  end  points  and  therefore 
the  terms  in  curly  brackets  must  be  identically  zero  yields 
an  equation  for  the  weighting  function  iji+(r,E,s): 


H+*  + (r ,E , s) 


^ »|/+(r,E,s) 


(5.6) 


and  requires  that  $(r,E,s)  satisfy  Eq.  (5.3).  The  stationary 
value  of  s is  given  by  the  following  expression  obtained 
from  Eq.  (5.4)  : 

/ dE  / dr  ij/+  (r,E,s  ) H !];  ( r rE , s ) 


/ dE  / dr  v (r ,E  ,sQ)  i *(r,E,so) 


(5.7) 


Equation  (5.7)  may  be  written  in  terms  of  reactivity,  mean 
generation  time,  and  effective  delayed  neutron  fraction 
as  follows: 
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P(s0)  - Beff(s0) 

A(so> 


(5.8) 


where  o(so),  6eff(sQ)  an<3  A(sq)  are  given  by 


/ dE  / dr  ij>+(r,E,s  ) L<ji(r,E,s  ) 

(so)  =1 2 2 

/ dE  / dr  ^+(?,E,so)x(E)P^(r/E' ,sQ) 


(5.9) 


Seff(so> 


D 

I Bi  / dE  / dr  (r,E,so)  fi(E)P4i(r,E',so) 
i=l 


/ dE  / dr  (r,E,sQ) x (E) P?(r,E' ,sq) 

(5.10) 


/ dE  / dr  £+(r,E,sQ)  i *(r,E,so) 

A(so)  » 

/ dE  / dr  * + (r,E,s0) x (E)P*(r,E* ,sq)  (5;ll) 

P*(r,E'(9o)  = / dE'  vEf(r,E'  ) ij»  ( r ,E',  sq) 

E'  (5.12) 

L?(r,E#sQ)  = - V*D(r,E)  7iF  (r,E,sQ)  + ET  (r  ,E)  (r  ,E  f sQ) 

/ dE'  E_(r,E'-»>E)*(r,E'  ,s  ) 

3 U 

E' 


(5.13) 
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Equation  (5.8)  follows  directly  from  Eq.  (5.7)  by 
rearranging  Eq.  (5,2)  in  terms  of  the  quantities  defined 
by.  Eqs.  ( 5 . 9)  - ( 5 . 13 ) . 

In  the  usual  1-point-model  formulation  of  Eq.  (5.9) 

separability  is  assumed  and  the  resulting  expressions  for 

P#  ®eff*  an<^  ^ are  ^dependent  of  s.  In  the  expressions 

for  p,  6 and  A the  functions  ^+(r,E,s  ) and  I(r,E.s  ) 

e11-  o o 

are  replaced  by  <J/(r,E,0)  which  satisfy  the  following 
equations 

Ho<Mr,E,0)  = 0 (5.14) 

H^+(£,E,0)  = o (5.15) 

where 

Ho<Mr,E,0)  = x(E)  / dE'v£f (r,E)*(r,E,0)  + 

E* 

7»D(r,E)  V<|>  (r,E,0)  - Et  (r,E)  <j>  (r,E  ,0)  + 

/ dE*  Eg  ( r /E ' -*-E)  ip  (r  fE'  , 0)  (5.16) 

E'  . 

The  difference  between  H as  given  by  Eq.  (5.2)  and  Hq  is 
in  the  fission  spectrum.  The  spectrum  of  all  fission  neutrons 
is  contained  in  Hq  whereas  only  the  prompt  neutron  spectrum 
reduced  by  a factor  of  (1-3)  is  contained  in  H. 
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In  the  numerical  computations  to  be  presented  later 
it  is  the  value  of  so  given  by  Eq.  (5.8)  that  is  referred  to 
as  the  "exact"  value  of  the  prompt  neutron  decay  constant. 


Computer  Codes 


The  objectives  for  performing  computations  on  the 
UFTR  were  the  following: 

(a)  to  obtain  a set  of  group  constants  that  well 
represent  the  neutronic  properties  of  the 
UFTR 

(b)  to  perform  a detailed  computation  of  the 
space-  and  energy-dependent  flux  and  adjoint 
flux 

(c)  to  compute  the  mean  generation  time  and  the 
effective  delayed  neutron  fraction 

(d)  to  compute  the  decay  constant  of  the  prompt 
neutron  die-away  using  the  exact  method  pre- 
sented in  the  preceding  section,  the  1-point 
approximation  and  several  coupled  approxima- 
tions and  to  compare  the  results  of  these 
computations  with  each  other  and  with  experi- 
mental measurements 


To  accomplish  these  objectives  it  was  necessary  to  utilize 

a series  of  computer  codes.  Figure  2 shows  a block  diagram 

of  the  codes  required  to  go  from  the  initial  basic  data 

to  the  final  result,  namely  the  decay  constant  for  the  prompt 

neutron  die-away.  One  additional  code  that  was  used  in  some 

» 

preliminary  investigations  is  the  TWENTY  GRAND.  The  compu- 
tations performed  by  these  codes  are  discussed  in  this  section. 


Group  Constants  — Two  codes,  TEMPEST  (Shu  60)  and 
GAM-I  (Joa  61) , were  used  to  compute  group  constants  for  the 
UFTR. 
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Theoretical  Computations 
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TEMPEST  is  a neutron  thermalization  code  based  upon 
the  Wigner-Wilkins  approximation  for  light  moderators  and 
the  Wilkins  approximation  for  heavy  moderators.  For  the 
UFTR  computations,  the  Wigner-Wilkins  approximation  was 
used  in  the  fuel  regions,  and  the  Wilkins  approximation 
was  used  in  the  graphite  regions. 

GAM-I  is  a consistent  multigroup  code  for  the 
calculation  of  fast  neutron  spectra  and  multigroup  constants. 
A very  complete  description  of  this  code  is  given  in 
reference  (Joa  61) . 

Flux  Computations  — The  principle  code  used  in  the 
flux  computations  was  the  SAIM  code.  SAIM  is  a modification 
of  the  FAIM  code  (Bal  62) , a multigroup,  1-dimensional  dif- 
fusion equation  code.  This  code  will  handle  as  many  as  18 
energy  groups,  101  space  points,  and  20  regions.  Neutrons 
are  permitted  to  scatter  from  a given  group  to  as  many  as 
8 lower  energy  groups.  Basically,  the  code  solves  the 
equation  given  by  Eq.  (3.1).  In  order  to  perform  the  compu- 
tations required  in  this  study  and  to  expedite  the  data 
processing  the  following  modifications  were  made  in  the 
FAIM  code : 

(1)  The  nuclear  data  library  was  removed 
because  it  was  not  necessary  and  only 
increased  the  running  time. 

, (2)  The  FAIM  code  is  a chain  job  and  there- 

fore cannot  be  run  under  monitor  control 
on  the  University  of  Florida  computer. 

After  removing  the  nuclear  data  library 
and  the  subroutines  that  used  this  data, 


\ 
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it  was  possible  to  reduce  the  size  of 
the  program  so  that  it  could  be  run 
under  monitor  control. 

(3)  The  FAIM  code  has  an  option  in  which 
it  will  search  for  the  transverse 
buckling  in  a region,  or  in  several 
regions,  required  to  make  the  reactor 
have  a specified  multiplication  constant. 
Thus,  the  FAIM  code  will  search  for  a 
region-dependent  transverse  buckling, 
but  not  a group-dependent  transverse 
buckling.  It  became  necessary  to  modify 
the  code  so  that  a group-dependent  trans- 
verse buckling  could  be  determined. 

(4)  Another  modification  was  the  option  to 
perform  the  partial  flux  and  partial  adjoint 
flux  computations  discussed  in  Chapter  III. 

(5)  To  expedite  data  processing  an  option 
was  added  for  punched  output  of  the  flux, 
adjoint  flux,  and  input  data.  This  punched 
output  may  be  read  directly  into  the 
C0UPKIN  code. 


Another  code  that  was  used  for  a preliminary  calcu- 
lation of  the  flux  is  the  TWENTY  GRAND,  a few-group,  2- 
dimensional  diffusion  equation  code  (Tob  62) . This  code 
is  capable  of  solving  neutron  diffusion  problems  in  cylin- 
drical or  slab  geometry  for  1 to  6 groups.  Up  to  3000  mesh 
points  may  be  used.  Execution  of  the  code  requires  an 
enormous  amount  of  computer  time.  For  example,  a 4 -group 
computation  with  a 50  by  50  mesh  point  layout  requires 
several  hours  of  running  time  on  the  IBM-709  computer. 
Unfortunately,  in  its  original  form  this  code  had  no  restart 
option  and  only  wrote  the  output  tape  at  the  completion  of 
the  computations.  As  a result  the  program  might  run  for 
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several  hours  without  any  indication  of  the  convergence  of 
the  calculation.  To  make  matters  worse,  if  the  desired 
convergence  criteria  were  not  met  before  the  maximum  number 
of  iterations  were  performed,  it  was  necessary  to  repeat 
all  of  these  iterations  if  it  were  later  decided  to  attain 
better  convergence  by  performing  more  interations.  Several 
modifications  were  made  in  the  code  to  remove  these  undesirable 
limitations.  The  modifications  are  summarized  as  follows: 

(1)  On-line  print-out  of  convergence 
information  occurs  every  time  the 
program  goes  through  L00P  iterations, 
where  L00P  is  a value  set  by  the  user. 

(2)  The  flux  is  stored  on  tape  at  the 
completion  of  each  L00P-th  iteration. 

(3)  The  user  has  the  option  to  obtain 
punched  output  of  the  flux  and 
adjoint  flux. 

(4)  A problem  may  be  restarted  either  with 
the  punched  output  or  with  the  storage 
tape. 

In  addition  to  group  constants,  flux,  and  a joint 

flux  it  is  also  necessary  to  compute  the  average  inverse 

speed  of  the  neutrons  in  each  group.  This  computation  re- 
% 

quired  another  computer  code  which  is  discussed  in  the 
next  section. 

Computation  of  the  Average  Inverse  Speed  — In  order 
to  compute  kinetics  parameters  it  is  necessary  to  determine 
the  average  inverse  speed  of  each  group.  A simple  example 
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will  serve  to  illustrate  why  the  average  inverse  speed  is 
needed  instead  of  the  average  speed.  The  mean  generation 
time  is  given  by 

/ dE  / dr  *+(?,E)  ^ y 

A = — - - 

/ dE  / dr  X(E)*+(?,E)  / dE'  vEf(r,E' )t(r,E' ) 

(5.17) 


Note  that  it  is  the  reciprocal  of  v(E)  that  appears  under 

« 

the  integral  in  the  numerator  of  Eq.  (5.17).  In  terms  of 
raultigroup  diffusion  theory  the  mean  generation  time  is 
given  by 


l ! a?  *+i(?)  (iy  ♦h?) 


I [ / d?  x* i*ti<?>  vihj) 

j=l  i=l 


(5.18) 


where 


1 ( x)  = the  flux  in  group  i at  point  r 

4>  + 1(r)  = the  adjoint^of  the  flux  in  group 

i at  point  r 


the  average  inverse  speed  of  the 
neutrons  in  group  i.  This  average 
is  also  space-dependent  if  the 
neutron  spectrum  is  space-dependent 
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Equation  (5.18)  indicates  that  the  generation  time  is 

directly  proportional  to  the  inverse  speed  of  the  neutrons. 

\ 

Therefore,  one  should  be  careful  in  his  selection  of  a 
method  of  averaging,  as  well  as  in  the  arithmetic  compu- 
tations. Consider  the  difference  between  / — \ and  -- 

\v/  <~v> 

for  neutrons  with  a Maxwell-Boltzmann  distribution, 
is  given  by 

1 


/ 


V (E) 


n (E)  dE 


/"  n (E)  dE 
o 


(5.19) 


and  <v>  is  given  by 


<,v> 


/"  v (E)  n (E)  dE 
o 

/"  n (E)  dE 
o 


(5.20) 


where 


n (E)  = M (E)  = (1) 


_L_  E*  e'EAT 


(kT)3/2 


neutrons  per  unit  energy  (Ros  61) 


(5.21) 


v (E)  = 1.38  x 109  y /IT 


cm/sec  (Gla  52) 


(5.22) 
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T = temperature,  °K 
k = the  Boltzmann  constant 
= 8.6275  x 10~5  ev/°K 

E = energy,  MeV 

Substituting  Eqs.  (5.21)  and  (5.22)  into  Eqs.  (5.19)  and 
(5.20)  and  performing  the  integration  yields: 


= i _L_ 

* <v) 


1.273 

<v> 


(5.23) 


It  is  seen  from  this  computation  that  a significant  error 

may  be  introduced  by  using  JL_  instead  of  / -).  The 

<v>  ' v/ 

result  of  this  example  justifies  the  somewhat  elaborate 
averaging  procedure  which  follows. 

Let  us  first  consider  the  average  inverse  speed  for 
the  thermal  group.  In  addition  to  the  thermal  group  con- 
stants the  TEMPEST  code  also  prints  out  values  for  the 
energy  flux  pi # where : 


the  flux  of  neutrons  per  unit  energy 
in  the  energy  interval  AE^ 


AE. 

X 


E.  - E,  , 

l i-1 


i > 1 


Making  the  substitutions  : 


n (E) 


♦ (E)  /v (E) 
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and 


v (E)  = 1.38  x 109  {e 


into  Eq.  (5.19)  yields: 


<0 


/ MAX  e”1^ (E)  dE 
o 


1.33  x 109  /'MAX  e-^  (e)  dE 


I (E"1)1  Fl  AEi 
i-1 


1.38  x 10y  N 


l (e^)1  f1  AE. 


i=l 


(5.24) 


EHAX  = the  uPPer  bound  of  the  thermal 

group 

I 

N = the  number  of  energy  increments 

into  which  the  range  from  zero 

to  is  divided 

MAX 


Maxwellian  averaged  values  are  used  for 

r-H\L 


E-1V  and 


<E‘7 


These  values  are  given  by  the  following  expressions: 


/ 1 M (E)  dE 


E 


-l\i 


Ei-1 


E . 

/ 1 M (E)  dE 


'i-1 


(5.25) 
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E 


/ 1 E_is  M (E)  dE 


<E'V>i 


Ei-i 


JEi  M(E)  dE 


Ei-1 


(5.26) 


Substituting  Eq.  (5.21)  into  Eqs.  (5.25)  and  (5.26)  and 
integrating,  yields 


(e-1)1  - 2 { 


E.  e 

l 


-E i /kT 


Ei-1  e 


-Ei_1A,r 


1 

+ — } 

kT 


( 


E 


-H\i 


kT  -E  . . /kT  -E  . /kT 

- 1 e 1*1  - e 1 } 


4i 


(5;27) 

(5.28) 


Ai  “ { 


-E./kT  -E  ,/kT 

Ei  6 1 + Ei-1  e 


} AEi 


(5.29) 


The  expressions  for  averaging  the  inverse  speed  of 
the  neutrons  above  thermal  are  in  terms  of  lethargy  rather 
than  energy.  In  addition  to  the  fast  group  constants  the 
GAM-I  code  also  prints  out  values  for  the  lethargy  flux. 


G , where 


the  flux  of  neutrons  per  unit 
lethargy  in  the  lethargy  interval 
AU . 


A°i 


u.  - u.  . 

1 1-1 


i > 1 
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The  relationship  between  lethargy  and  energy  is  given  by 

-U 

E = 10  e 


and  the  speed  of  a neutron  with  lethargy  U is  expressed  by 


V(U) 


/lcT  x 1.38  x 109 


-U/2 

e 


In  terms  of  lethargy  the  average  inverse  speed  for  some 
I-th  fast  group  is  given  by 


/,Ui  vW  n(u)  du 

/1\I  _ Ui~l 

\ v)  = " 

JUi  n (U)  dU 


'10  x 1.38  x 10' 


U.  U 

/ 1 e <fr(U)  dU 


U/2 

e <fr(U)  dU 


U 


i-1 


a 


i g1  aui 

iel 


1.38  x 109 


U/2V  . 

I (e  y G1  AUX 
iel 


(5.30) 
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Lethargy  averaged  values  are  used  for  (^eU^)  * and  <eU/2>i- 
These  values  are  given  by  the  following  expressions 

/ 1 eU  dU 


<euy 


Vi 


AU. 


U.  U.  , 
e 1 - e 1“1 


AU , 


(5.31) 


and  similarly, 


(eU/2)i 


U./2  Ui-1/2 

2 (e  1 " e ) 


AU. 

1 


(5.32) 


A code,  called  EV,  has  been  written  to  compute  values 
for  the  average  inverse  speed  using  the  method  discussed  in 
this  section.  In  addition  to  the  inverse  speed  the  program 
also  computes  the  average  energy  of  each  group. 


Computation  of  the  Kinetics  Parameters  — A code 
has  been  written  to  perform  the  computation  of  kinetics 
parameters  according  to  the  multigroup  formulation  presented 
in  Chapter  III.  The  name  of  the  code  is  C0UPKIN.  It  is 
designed  to  take  the  punched  output  from  the  SAIM  code 
without  additional  processing.  The  code  is  dimensioned  to 
compute  the  coupling  parameters  for  a reactor  described  by 
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up  to  20  broad  regions  and  up  to  18  broad  groups. 

It  will  perforin  computations  for  1-dimensional  slab, 
spherical  or  cylindrical  geometries.  Numerical  integra- 
tions are  performed  with  the  use  of  a modified  trapezoidal 
rule  given  on  page  118  of  reference  (Mil  49) , and  numerical 
derivatives  are  computed  with  a 5— point  scheme  given  on 
page  97  of  reference  (Mil  49) . 

Eigenvalue  Computations  — The  eigenvalues  in  this 
study  were  computed  with  a subroutine,  EIG4,  obtained  from 
Ray  Booth.  The  subroutine  was  originally  obtained  through 
SHARE  and  has  the  identification  name  NU  EIG4  (Par  63). 

This  subroutine  computes  the  eigenvalues  of  a given  complex 
h x N matrix  for  2<N<70.  It  was  only  necessary  for  the 
author  to  write  a calling  program  in  order  to  use  this  sub- 
routine . 

V 

Results  of  UFTR  Computations 

Figure  3 shows  a schematic  drawing  of  the  core  of 
the  UFTR.  A more  detailed  description  of  the  core  is  given 
in  Thompson's  thesis  (Tho  61) . All  of  the  computations  pre- 
sented in  this  study  are  based  on  a 3500  g U-235  loading 
of  20  per  cent  enriched  uranium.  The  atom  densities  for 
this  arrangement  are  given  in  Table  A-l  of  the  Appendix. 

For  this  loading  the  value  of  as  determined  by  positive 

period  measurements  is  in  the  range  of  1.0030  to  1.0057 
depending  on  the  fuel  arrangement,  water  temperature,  and 
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152.5  cm 


Two-Dimensional  Representation 
b c c b 
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m Graphite 
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a = 47.6  cm 
b = 13.4  cm 
c =»  15.2  cm 
d = 49.9  cm 
e = 13.2  cm 


Section  E-E 

One-Dimensional  Representation 
Of  UFTR 


Figure  3.  Schematic  Drawing  of  the  University 
of  Florida  Training  Reactor 


109 


any  extraneous  poison  that  might  be  present  in  the  core 
(But  64).  In  this  study  a value  of  1.0057  was  assumed  for 
with  all  control  blades  removed  and  with  a core 
temperature  of  300°I<. 

Before  undertaking  an  18-group  calculation,  4-group 
calculations  were  made  with  the  SAIM  and  TWENTY  GRAND  codes. 

The  4-group  SAIM  results  were  checked  against  the  output  of 
the  old  version  of  the  FAIM  code  for  verification  that  the 
SAIM  code  was  performing  properly.  The  4-group  TWENTY  GRAND 
output  was  usqd  to  obtain  an  estimate  of  the  buckling.  The 
original  group  constant  calculations  were  based  on  experi- 
mentally measured  values  of  the  buckling  reported  in  Thompson's 
thesis.  Then,  after  the  4-group  computations  were  made,  the 
estimated  values  of  the  buckling  were  revised.  On  the  basis 
of  hand  calculation  using  the  output  from  TWENTY  GRAND,  the 
values  given  in  Table  3 were  used  in  the  final  computation 
of  the  constants  for  18  groups.  In  computing  values  for  the 
total  buckling  of  a given  region  it  was  assumed  that  the 
buckling  in  the  vertical  direction  equals  that  in  the  Y 
direction.  The  X and  Y directions  are  shown  in  Figure  3. 

The  vertical  direction  is  normal  to  the  X-Y  plane.  The 
group  constants  used  in  the  4-group  computations  are  listed 
in  Butterfield's  thesis  (But  65).  For  the  18-group  compu- 
tations the  constants  are  listed  in  Tables  A-2  through  A-6 
in  the  Appendix.  The  values  of  x*"  given  in  Table  A-5 
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TABLE  3 

BUCKLING  VALUES  USED  IN  GROUP 
CONSTANT  CALCULATIONS 


Group  

Regions  1 & 

2 -2 
a (cm  *) 

5 Regions  2 & 4 Region  3 

1,2, 3,4  -0.0101 

0.0112  -0.00304 

5,6,7,8,9,10,11  -0.00516 

0.00421  0.00239 

12,13,14,15,16,17  -0.00145 

0.00312  0.00250 

18  0.00282 

0.00222  0.00282 

Ill 


were  obtained  from  values  taken  from  ANL-5800  (Arg  63) 
manipulated  in  the  following  manner. 

Table  1-11  of  ANL-5800  gives  the  U-235  prompt 
fission  spectrum  as  a function  of  lethargy  and  Table  1-17 
gives  the  mean  energies  of  the  delayed  groups.  In  terms  of 
the  prompt  fission  spectrum  and  the  delayed  fission  spectrum 

the  total  fission  spectrum  is  given  by 

6 

X(U)  = (1-8)  F (U)  + l 6(U-U..)8. 

j = l 

where 

F(u)  = the  prompt  fission  spectrum 

0 . = the  average  lethargy  of  the 

-1  j-th  delayed  group 

The  value  of  x1  in  Table  A-5  is  given  by 


i 

X 


U. 

X1  - / 1 X (U)  du 

Ui-1 


u.  D u. 

(1-8)  / 1 F (U)  dU  + l 8^  / 5(U-U..)  dU 


U 


i-1 


j = l 


U 


i-1 


(5.33) 


The  integral  over  F(U)  is  given  in  Table  1-11  and  mean 
energies  corresponding  to  0^  are  given  in  Table  1-17. 

The  EV  code  was  used  to  compute  the  average  inverse 
speed  and  average  energy  of  the  neutrons  in  each  group  and 
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these  values  are  given  in  Tables  A-7  and  A-8  of  the 
Appendix.  Since  the  neutrons  in  the  thermal  group  have  a 
Maxwell-Boltzmann  distribution,  we  may  use  the  following 
formula  (Ros  61)  based  on  the  kinetic  theory  of  gases  to 
compute  an  effective  neutron  temperature  (Meg  60) : 

(E)  = \ k Tn 

where  (e^  is  given  in  Table  A-8  and  k is  the  Boltzmann 
constant.  The  effective  neutron  temperature  for  the  thermal 
neutrons  in  the  fuel  region  is  153°C  and  in  the  graphite  it 
is  93°C.  On  the  bhsis  of  differential  measurements  of  the 
spectrum  at  the  center  of  the  UFTR  Salah  determined  a 
value  of  115  * 6°C  for  T^  and  using  an  integral  method 
he  found  a value  of  76  * 2.9°C  (Sal  64).  The  theoretical 
value  of  93°C  is  bracketed  by  the  two  measured  values. 

In  performing  a 1-dimensional  calculation  of  a 3- 
dimensional  core  it  is  necessary  to  account  for  neutron 
diffusion  in  the  transverse  direction  by  a buckling  term. 

If  the  flux  in  the  transverse  direction  in  known,  the 
buckling  can  be  computed  according  to  the  formula 

/ dE  / dr  <j>f(?,E)  V2  * (r,E) 
j y 

B = 

/ dE  / dr  $+(r,E)  $(r,E) 

where  V2  is  the  Laplacian  in  the  Y-Z  plane.  However,  in 
Y • z 

2 

these  calculations  since  the  value  of  B was  not  known  very 
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accurately  as  compared  to  the  value  of  k a search  was 

made  by  the  SAIM  code  for  the  value  of  the  transverse  buckling 

that  would  yield  a value  of  1.0057  for  k A constant 

eft 

value  was  assumed  for  all  groups  and  all  regions.  The 

. —2  ' 
search  yielded  a value  of  0.00254  cm  for  the  transverse 

buckling.  A search  was  then  made  to  determine  the  amount  of 

thermal  group  poison  that  must  be  added  to  region  2 in  order 

to  reduce  k^^  from  1.0057  to  1.0000.  This  search  yielded 

a value  of  0.00109  cm  ^ for  the  macroscopic  cross-section 

of  the  poison.  The  C0UPKIN  code  was  then  used  for  computing 

3eff  an^  A averaged  over  all  groups  and  all  space  points. 

The  resulting  values  were  8 = 0.00696  and  A = 188  ysec 

Gtt 

which  gives  a value  of  37.0  sec-1  for  8 

ef  f 

According  to  the  1-point-model  approximation,  the 
decay  constant  of  the  prompt  neutron  die-away  for  a reactor 
pulsed  at  delayed  critical  is  given  by  the  ratio  of  8^^/A. 
However,  according  to  Butterfield's  measurements  (But  65) 
the  decay  constant  at  delayed  critical  is  22.2  * 0.5  sec”1, 
and  according  to  random  noise  measurements  made  by 
Kavipurapu  (Kav  64) , the  value  is  29  * 5 sec-1  (the  error 
of  *5  was  suggested  by  Ohanian  (Oha  65)  because  of  poor 
statistics  in  the  data) . It  was  at  first  thought  that  the 
wide  discrepancy  between  the  theoretical  value  and  the 
experimental  value  was  due  to  the  inability  of  1-point-model 
kinetics  to  accurately  predict  the  decay  constant  for  a 
reactor  with  weakly  coupled  fuel  regions.  However,  when 
an  exact  computation  of  the  decay  constant  was  found  to  be 
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within  a few  per  cent  of  the  1-point-model  value,  it  became 
apparent  that  the  discrepancy  must  be  due  to  some  computa- 
tional error. 

The  simple  problem  that  was  analyzed  in  Chapter  IV 
gave  a clue  as  to  where  a major  portion  of  the  apparent 
error  might  result.  Consider  the  thermal  group.  The  mean 
generation  time  averaged  over  all  groups  and  all  regions  is 
strongly  dependent  on  the  weighted  neutron  population  in 
the  thermal  group.  Since  most  of  the  thermal  neutrons  are 
in  the  graphite,  then  the  mean  generation  time  is  strongly 
dependent  on  the  removal  cross-section  in  the  thermal  group. 
The  thermal  removal  cross-section  is  given  by 


For  the  graphite  E^h  = 0.000241  cm-1,  Dth  = 0.9156  cm  and 

a 

if  = 0.00254  cm  then 

- 0.000241  + 0.00233  = 0.00257  cm"1 


In  this  case  approximately  90  per  cent  of  the  thermal  re- 
moval in  the  graphite  is  attributed  to  transverse  leakage. 

Let  us  examine  Figures  4 through  7.  These  figures 
indicate  that  the  distribution  of  the  flux  is  not  strongly 
region-dependent;  however,  observe  the  change  in  curvature 
from  one  group  to  the  next.  To  emphasize  the  difference 
in  curvature  and  to  obtain  an  approximate  value  for  the 
transverse  buckling,  each  curve  has  been  fitted  with  a 


Flux 
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Figure  4 . Group  1 Flux  as  Computed  by 
the  TWENTY  GRAND  Code 


Flux 
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Figure  5.  Group  2 Flux  as  Computed  by 
the  TWENTY  GRAND  Code 


Flux 
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Figure  6.  Group  3 Flux  as  Computed  by 
the  TWENTY  GRAND  Code 


Flux 
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Figure  7.  Group  4 (Thermal)  Flux  as  Computed 
by  the  TWENTY  GRAND  Code 
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cosine  shape.  Values  were  obtained  for  the  transverse 
buckling  for  the  top  throe  groups  by  assuming  that  the 
flux  has  the  same  shape  in  the  vertical  direction  as  in 
the  Y direction.  For  the  thermal  group  the  actual  height 
of  the  graphite,  142  cm,  was  used  to  compute  the  vertical 
buckling  which  was  added  to  the  buckling  in  the  Y direction 
to  obtain  the  transverse  buckling.  These  values  are  listed 
in  Table  4 . 

Another  flux  computation  was  performed  with  the  SAIM 

code.  In  this  computation  the  transverse  buckling  for  the 

_2 

thermal  group  was  held  constant  at  0.000921  cm  and  a 

search  was  made  in  groups  1 through  17  that  would  yield  a 

-2 

value  of  1.0057  for  k ff.  A value  of  0.00304  cm  was 
determined  for  the  transverse  buckling.  A search  was  then 
made  to  determine  the  amount  of  thermal  group  poison  that 
must  be  added  to  region  2 in  order  to  reduce  k^^  from 
1.0057  to  1.0000.  This  search  yielded  a value  of  0.000909  cm 
for  the  macroscopic  cross-section  of  the  poison.  When  the 
C0UPKIN  code  was  used  to  compute  and  A,  their  values 

were  determined  to  be  $eff  ~ 0.00703  and  A =>  247.5  ysec 
which  gives  a value  of  28.4  sec  ^ for  Be^f/A.  This  value 
agrees  reasonably  well  with  the  experimentally  measured 
values  given  earlier.  Also,  Chastain  (Cha  61)  deduced  from 
transfer  function  measurements  a value  of  240  ysec  for  the 
prompt  neutron  generation  time  which  is  in  good  agreement 
with  the  calculated  value  of  A. 
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TABLE  4 


GROUP-DEPENDENT  VALUES  OF  THE  TRANSVERSE 
BUCKLING  OBTAINED  FROM  THE  TWENTY 
GRAND  COMPUTATIONS 


Group 

Transverse  3uckling 

(Based  on  a 4-qroup 
Calculation) 

, -2, 
(cm  ) 

1 

0.00309 

2 

0.00245 

3 

0.00182 

4 

0.000921 
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. On  the  basis  of  the  preceding  comparison  of 
calculated  values  with  measured  values,  it  appears  that 
the  computed  group  constants  give  a fairly  accurate  de- 
scription of  the  neutronic  properties  of  the  UFTR  in  terms 
of  multigroup  diffusion  theory.  Certainly,  the  breach 
between  theory  and  experiment  has  been  closed,  e.g.,  the 
prompt  neutron  generation  time  given  in  the  hazards  report 
for  the  UFTR  (Dun  58)  is  140  ysec  which  is  far  below  any 
values  deduced  from  experimental  measurements.  Of  course, 
refinements  can  be  made  in  these  flux  calculations  and 
several  are  suggested  later  in  this  section. 

To  check  the  accuracy  of  the  1-point  model  an  exact 
computation  of  the  decay  constant  was  performed.  Also, 
several  coupled  representations  were  included  in  the  study. 
The  results  of  these  computations  are  listed  in  Table  5. 
Values  for  the  coupled  kinetics  parameters  are  given  in 
Tables  6,  7 and  8,  The  nomenclature  for  the  coupled  repre- 
sentation is  as  follows:  for  the  2— broad— region  representa- 

tion, broad-region < 1 extends  from  A to  D shown  in  Section 
E-E  of  Figure  3 and  broad-region  2 extends  from  D to  G; 
for  the  3-broad-region  representation,  broad-region  1 extends 
from  \ to  C,  broad-region  2 extends  from  C to  E and  broad- 
region  3 extends  from  E to  G;  and  for  the  2-broad-group 
representation,  broad-group  1 includes  the  first  seventeen 
grouos  and  broad-grouo  2 is  the  thermal  group.  The  use  of 
these  parameters  for  studying  space-dependent  kinetics  is 
discussed  in  a later  section. 


TABLE  5 


VALUES  OF  THE  DECAY  CONSTANT  FOR  THE 
PROMPT  NEUTRON  DIE-AWAY 


Model 

Number  of 
Points  in 
Phase  Space 

Decay  Constant 

, -In 

( sec  ) 

Exact 

_ 

-27.6 

1-broad-region , 

1-broad -group 

1 

1 

to 

00 

• 

2-broad-regions , 

1-broad-group 

2 

-27.6 

3-broad-regions , 

1-broad-group 

3 

-28.4 

2-broad-regions , 

2 -broad-groups 

4 

o 

• 

CO 

CN| 

1 

3-broad-regions , 

2 -broad-groups 

6 

-27.9 

KINETICS  PARAMETERS  FOR  A 2-REGION  AN! 
3-REGION  REPRESENTATION  OF  THE  UFTR 
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Note:  Number  in  parentheses  indicates  the  power  of  10  by  which  the  given 

value  is  to  be  multiplied. 


KINETICS  PARAMETERS  FOR  A 2 -REGION 
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In  Table  5 the  values  given  by  the  coupled  point-model 

representation  are  only  good  to  two  significant  figures.  The 

determinants  from  which  these  eigenvalues  are  computed  have 

3 6 

elements  of  the  order  of  10  to  10  and  although  each  element 
is  accurate  to  four  significant  figures,  accuracy  is  lost 
in  taking  small  differences  between  large  numbers.  The 
next  to  the  smallest  eigenvalues  in  each  of  the  coupled 
representations  is  approximately  -5000  sec”''".  In  all  cases 
all  of  the  eigenvalues  are  pure  real.  The  numerical  accu- 
racy in  determining  the  decay  constant  from  the  exact  compu- 
tation and  the  1-point-model  computation  is  much  better 
because  both  of  these  methods  only  require  taking  the  ratio 
of  two  values  which  are  accurate  to  four  significant  figures. 
In  comparing  the  value  given  by  the  exact  computation  with 
that  given  by  the  1-point  model,  it  is  observed  that  they 
differ  by  only  2.9  per  cent. 

The  agreement  between  the  exact  computation  and  the 
1-point  approximation  as  well  as  the  coupled  approximations 
is  very  good.  Much  of  the  difference  between  the  theoretical 
computation  and  the  experimental  measurements  can  be  attri- 
buted to  the  difficulty  in  representing  the  UFTR  by  a 1- 
dimensional  model.  There  are  some  definite  improvements 
that  can  be  made  in  the  theoretical  computations.  The 
following  suggestions  are  offered  to  anyone  who  might  find 
it  necessary  to  refine  the  computations. 
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(1)  The  TWENTY  GRAND  calculations  should  be 
repeated  using  a refined  set  of  cross- 
sections  and  group-dependent  values  for 
the  vertical  buckling.  Table  9 gives 
the  recommended  energy  boundaries  for  a 
4-group  calculation/  the  total  buckling 
to  be  used  in  the  group  constant  calcu- 

• lations,  and  the  group-dependent  vertical 
buckling  for  the  TWENTY  GRAND  calculation. 

For  these  computations  symmetry  may  be 
assumed  and  only  one  quadrant  of  the  core 
need  be  considered.  This  TWENTY  GRAND 
calculation  would  require  about  two  hours 
of  time  on  the  University  of  Florida  IBM- 
709  computer  or  about  one  hour  on  the 
CDC-1604  computer  at  the  Oak  Ridge  National 
Laboratory.  If  the  problem  is  run  on  the 
CDC-1604  computer,  punched  output  should 
be  requested. 

(2)  The  revised  TWENTY  GRAND  calculations  may  be 
used  to  bbtain  better  values  of  the  transverse 
buckling  and  the  total  buckling.  These  re- 
vised values  for  the  buckling  based  on  the 
4-group  calculation,  may  then  be  used  to 
compute  a new  set  of  group  constants  for 

an  18-group  calculation. 

(3)  The  values  given  in  Table  7 show  that  97.3 
per  cent  of  the  weighted  neutron  population, 
i.e,,  $+ 1/v  <j> , is  in  the  thermal  group.  Be- 
cause of  the  significance  of  the  thermal 
neutrons  in  the  UFTR,  one  should  give  some 
consideration  to  performing  a multigroup 
calculation  with  several  thermal  groups. 


In  the  course  of  this  study  several  volumes  of 
numerical  data  have  accumulated.  The  computer  output  has 
been  bound  and  labeled  and  is  available  to  anyone  who  re- 
quests it.  Also,  the  group  constants,  the  output  from  the 
TWENTY  GRAND  calculations,  and  the  output  from  the  SAIM 
runs  are  on  punched  cards.  Given  the  format,  a user  may 
process  this  information  to  suit  his  needs.  The  output  from 


RECOMMENDED  BUCKLING  VALUES  FOR 
FUTURE  CALCULATIONS 
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the  final  SAIM  run  is  summarized  in  Figures  8 through  11. 

The  values  in  these  graphs  were  obtained  using  a transverse 

-2 

buckling  of  0.000921  cm  for  the  thermal  group  and  0.00304 

_ 2 

cm  for  the  other  groups  with  a thermal  poison  cross-section 
of  0.000909  cm”1  in  region  2.  The  values  in  these  four 
figures  all  have  the  same  normalization. 

The  results  of  the  numerical  computations  presented 
in  this  section  have  shown  that  the  1-point  model  is 
adequate  to  predict  the  asymptotic  decay  constant  measured 
in  a pulsed  neutron-beam  experiment.  Therefore,  it  is  not 
necessary  to  use  a more  elaborate  kinetics  model  to  analyze 
this  experiment.  However,  if  it  is  the  transfer  function 
that  is  being  measured,  then  it  is  necessary  to  use  a set 
of  coupled  point  models  whose  behavior  is  described  by  the 
set  of  coupled  kinetics  equations  derived  in  Chapter  II. 

In  the  next  section  it  is  shown  how  the  numerical  values 
given  in  Tables  6,  7 and  8 may  be  used  to  predict  the 
response  of  one  region  of  the  reactor  due  to  a time-dependent 
perturbation  in  another  region. 

Parameters  for  a Space-Dependent  Kinetics  Study 

In  the  section  entitled  Application  to  the  Pile 
Oscillator  Experiment  in  Chapter  II  a region-  and  group- 
dependent  transfer  function  was  derived.  Values  for  all 
of  the  parameters  appearing  in  this  transfer  function  are 
given  in  Tables  6,  7 and  8 for  a 1-mode  representation  of 
the  flux.  If  the  input  perturbation  is  a time  varying 


Flux  or  Adjoint  Flux 
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Figure  8.  Flux  and  Adjoint  Flux  for  the 

Thermal  Group  as  Computed  by  SAIM 


Flux  or  Adjoint  Flux 
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X (cm) 


Figure  9.  Flux  and  Adjoint  Flux  Summed 
Over  the  Groups  Above  Thermal 
as  Computed  by  SAIM 


(x)  <p  (x) 
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X (cm) 


Figure  10.  Adjoint  Weighted  Flux  for  the 
Thermal  Group 


(x)  <t>  (x) 
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Figure  11 


Adjoint  Weighted  Flux  for  the 
Neutrons  Above  Thermal 
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thermal  group  poison,  the  parameter  c^,,(t)  representing  the 
input  perturbation  is  given  by:  y 


ek-(t> 

y' 


{/  dE  / dr  4 + $E  }f(t) 

E^,  Vol.  of  poison  p 

{(DV4  + .V* *)k,  + (*  + ET4>)k,  + y / dE  / dsk, 


j / dE  / dsk,  O^D74>*nk  } 


(5.34) 


In  Eq.  (5.34)  y'  represents  the  thermal  group  if  a 2-broad- 
group  representation  is  used.  The  subscript  k'  represents 
the  broad  region  in  which  the  perturbation  is  located.  Tfye 
function  f(t)  gives  the  time-dependence  of  the  poison  with 
a magnitude  of  Ep.  Values  for  the  denominator  of  the 
expression  on  the  right-hand  side  of  Eq.  (5.34)  are  given 

• f 

m Table  10  and  values  of  $ $ for  the  thermal  group  are 
given  in  Figure  10.  In  this  1-dimensional  representation 
of  the  UFTR  it  is  not  possible  to  calculate  E D directly 
but  it  can  be  determined  indirectly  by  a combination  of 
reactivity  measurements  and  eigenvalue  calculations.  With 
the  poison  in  a given  position  bring  the  reactor  to  critical 
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TABLE  10 


DENOMINATOR  OF  THE  EXPRESSION  ON  THE 
RIGHT-HAND  SIDE  OF  EQUATION  (5.34) 


Model 

Region 

Denominator 

2-broad-regions , 

1-broad-group 

1 

5.5560737 

2 

6.5642879 

2-broad  regions, 

2 -broad-groups 

1 

0.6864955 

2 

0.7825920 

3-broad-regions , 

1-broad-group 

1 

4.8908975 

2 

2.3391127 

3 

5.9057559 

3-broad-regions , 

2 -broad-groups 

1 

0.6320942 

2 

0.3863222 

3 

0.7489065 
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and  determine  the  worth  of  the  poison  in  that  position  by 

referring  to  a control  rod  calibration  curve.  Then, 

perforin  a SAIM  calculation  and  have  the  program  search 

for  the  value  of  the  thermal  group  poison  in  the  given 

position  that  will  produce  the  measured  reactivity.  By 

repeating  this  procedure  for  several  positions  of  the 

poison,  it  will  be  possible  to  determine  E as  a function 

P 

of  position.  The  function  f(t)  gives  the  position  as  a 

function  of  time,  so  the  product  of  f(t)  and  l gives  the 

P 

variation  of  the  poison  as  a function  of  time. 

Coupled  kinetics  parameters  for  a 2-mode  representa- 
tion of  the  flux  will  not  be  given.  Early  in  the  study  when  it 
was  thought  that  the  discrepancy  between  the  theoretical 
computations  and  the  experimental  measurements  of  the 
asymptotic  decay  constant  was  due  to  the  weak  coupling 
between  the  fuel  regions,  a 2-mode  calculation  of  the  flux 
was  performed.  Figure  12  shows  a plot  of  the  total  flux 
and  the  partial  fluxes  resulting  from  an  18— group  computa- 
tion by  the  SAIM  code.  When  the  asumptotic  decay  constant 
was  determined  for  several  2-mode,  multi-broad-region , multi- 
broad-group  representations,  the  values  that  resulted  were 
in  good  agreement  with  the  1-point  value.  After  the  trans- 
verse buckling  values  were  revised,  the  2-mode  calculations 
were  not  repeated.  These  computations  are  expensive  and  it 
was  felt  that  the  simpler  1-mode  representation,  multi-broad- 
region,  multi-broad-group  representation  should  be  studied 


Flux 
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X (cm) 


Figure  12.  Flux  in  the  UFTR  as  Computed  by  the  SAIM  Code; 

is  the  Flux  Due  to  Neutrons  Born  in  the 

i-th  Region 
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more  thoroughly  before  going  to  a more  elaborate  model. 
However , it  should  be  noted  that  the  codes  for  performing 
the  modal  calculations,  namely  SAIM  and  C0UPKIN,  are 
already  written. 

Conclusions 

The  results  of  computer  calculations  of  the  flux 
in  the  University  of  Florida  Training  Reactor  have  been 
presented  along  with  a brief  discussion  of  the  codes  used 
to  make  these  computations.  Theoretical  computations  of 
the  asymptotic  decay  constant  of  the  prompt  neutron  die- 
away  at  delayed  critical  agreed  favorably  with  the  values 
determined  by  a pulsed  neutron-beam  experiment  and  a random 
noise  experiment.  Furthermore,  it  was  shown  that  the  1- 
point-model  kinetics  equation  very  accurately  predicts  the 
value  of  this  decay  constant.  For  the  UFTR  the  values 

0.007-03  and  A=  247.5  psec  were  computed.  A procedure 
for  improving  these  values  was  suggested. 

In  addition  to  the  results  of  the  flux  calculations, 
kinetics  parameters  have  been  presented  for  several  different 
coupled  point-model  representations  of  the  UFTR  and  a pro- 
cedure for  using  these  parameters  to  determine  a space- 
dependent  transfer  function  has  been  suggested. 


CHAPTER  VI 


SUMMARY  AND  CONCLUSIONS 

A general  method  of  studying  energy-  and  space- 
dependent  reactor  kinetics  problems  in  terms  of  a nodal,  or 
point-model,  representation  has  been  presented  and  estab- 
lished on  a firm  basis  by  deriving  the  pertinent  set  of 
equations  from  the  general  integro-dif ferential  form  of 
the  Boltzmann  neutron  transport  equation.  The  formulation 
of  the  equations  from  first  principles  on  the  basis  of  a 
semidirect  variational  method  guided  the  defining  of 
coupling  parameters  in  a consistent  way.  The  derivation 
was  so  general  as  to  permit  the  reduction  of  the  space-  , 
energy-,  angle-  and  time-dependent  reactor  neutron  field 
to  any  number  of  nodes  that  depend  only  on  time.  These 
nodes  may  represent  the  integrated  neutron  field  divided 
into  energy  groups,  configuration  regions,  or  a combination 
of  the  two.  Furthermore,  the  steady-state  flux  and  adjoint 
flux  from  which  the  coupling  parameters  are  computed  may  be 
represented  by  any  number  of  eigenmodes. 

For  the  sake  of  generality  the  coupled  kinetics 
equation  was  derived  in  terms  of  the  neutron  transport 
equation;  however,  in  order  to  facilitate  the  computation 
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of  the  coupling  parameters  and  to  clarify  the  general 
treatment,  formulas  were  given  for  computing  the  coupling 
parameters  in  terms  of  diffusion  theory. 

In  connection  with  this  work  a number  of  computer 
codes  were  developed  for  use  on  the  University  of  Florida 
IBM-709  computer.  These  codes  have  a wide  range  of  appli- 
cability for  studying  both  static  and  dynamic  reactor 
physics  problems  and  are  valuable  additions  to  the  research 
facilities  of  the  Nuclear  Engineering  Department. 

Both  a 4-group,  2-dimensional  and  an  18-group,  1- 
dimensional  flux  calculation  for  the  University  of  Florida 
Training  Reactor  have  been  performed.  A summary  of  these 
computations  was  presented  in  Chapter  V.  A complete  listing 
of  the  printed  output  has  tbeen  bound  under  separate  cover. 
Also,  punched  output  from  these  computer  runs  is  filed  in 
the  IBM  card  storage  room  of  the  Nuclear  Engineering  Depart- 
ment. The  cross-sections  for  the  18-group  calculations  are 
listed  in  the  Appendix  and  those  for  the  4-group  calculations 
are  listed  in  Butterfield's  thesis  (But  65).  This  infor- 
mation will  be  valuable  for  future  theoretical  studies  of 
the  UFTR. 

An  exact  computation  of  the  asymptotic  prompt  neutron 
decay  constant  was  used  to  check  the  accuracy  of  the  theory 
presented  in  Chapter  II.  The  decay  constant  predicted  by 
several  point-models  ranging  from  a 1-point  model  to  a 
6-point  model  (3  regions  and  2 groups)  was  computed.  The 
value  predicted  by  all  of  the  point  models  was  within 
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3 per  cent  of  the  value  predicted  by  the  exact  computation. 
The  results  of  this  comparison  give  credence  to  the  validity 
of  the  coupled  kinetics  formulation.  Another  important 
result  of  this  comparison  was  the  close  agreement  of  the 
values  predicted  by  the  1-point-model  computation  and  the 

i 

exact  computation.  This  result  indicates  that  the  1-point- 
model  kinetics  equation  is  sufficiently  accurate  to  analyze 
neutron  die-away  measurements  even  in  a reactor  with  weakly 
coupled  fuel  regions. 

The  following  values  were  computed  with  the  18- 

group,  1-dimensional  flux  and  adjoint  flux: 

6 =>  0.00703 

ef  r 

A =*  247.5  usee 

The  exact  computation  of  the  asymptotic  prompt  neutron  decay 

constant  yielded  a value  of  -27.6  sec  ^ which  is  in  fair 

agreement  with  measured  values  of  -22.2  * 0.5  (But  65)  and 

-29  * 5 sec”1'  (Kav  64)  . Suggestions  were  made  for  refining 

the  theoretical  computations.  Of  the  two  parameters,  Be££ 

and  A,  the  value  for  3 is  by  far  the  more  accurate.  It 

ef  f 

was  shown  in  Chapters  IV  and  V that  the  value  for  A is  a 
sensitive  function  of  the  transverse  buckling  which  is  not 
very  accurately  known  for  the  UFTR. 

Multiple  point-model  representations  of  the  reactor 
neutron  field  should  find  their  greatest  applicability  in  the 
study  of  transients  and  of  the  propagation  of  perturbations 
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in  a reactor  with  weakly  coupled  fuel  regions.  A general 
formula  was  derived  in  Chapter  II  for  a region-  and  group- 
dependent  transfer  function  and  values  for  the  parameters 
appearing  in  this  transfer  function  were  reported  in 
Chapter  V with  the  LJFTR  data. 


APPENDIX 


APPENDIX 


Presented  on  the  following  pages  are  tables 
describing  the  neutronic  properties  of  the  UFTR  for 
an  18-group  diffusion  theory  model.  These  values  are 
based  on  a 20  per  cent  enriched  uranium  fuel  loading 
with  a 3500  g U-235  critical  mass. 
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REGION 


2 AND 
2 AND 
2 AND 
2 AND 
li  3 AND 


NOTE 


TABLE  A-l 
ATOM  DENSITIES 
FOR  THE  UFTR  ( THO  61) 
WITH  20  PERCENT  0-235  FUEL 


ELEMENT  DENSITY 

(ATOMS/CC) 


HYDROGEN 
U-235 
U-238 
ALUMINUM 
5 GRAPHITE 


3.9141+22) 
1 . C69 ( +20 ) 
4.222 (+20) 
1 . 997 ( +22 ) 
8 . 02  5 ( +22 ) 


NUMBER  IN  PARENTHESES  INDICATES  THE  POWER  OF  10 
BY  WHICH  THE  GIVEN  VALUE  IS  TO  BE  MULTIPLIED. 
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TABLE  A- 2 
GROUP  BOUNDARIES 
BY  ENERGY  AND  BY  LETHARGY 


GROUP 

ENERGY  (MEV) 

LETHARGY 

NUMBER 

LOWER  BOUNDARY 

UPPER  BOUNDARY 

10.0000000 

0.0 

1 

6.0653066 

0.5 

2 

3.6787944 

1.0 

3 

2.2313016 

1.5 

4 

4.97870681-1) 

3.0 

5 

1.83156391-1) 

4.0 

6 

1.  11089961-1) 

4.5 

7 

6.73794701-2) 

5.0 

8 

4.08677151-2) 

5.5 

9 

2.47875221-2) 

6.0 

10 

1.50343921-2) 

6.5 

1 1 

9.11881971-4) 

7.0 

12 

5.53084371-4) 

7.5 

13 

3.35462631-4) 

8.0 

14 

1.23409801-4) 

9.0 

15 

4.53999301-5 ) 

10.0 

16 

6.14421231-6) 

12.0 

17 

6.82560351-7) 

16.5 

18 

0.0 

INFINITY 

NO  I E 


NUMBER  IN  PARENTHESES  INDICATES  THE  POWER  OF  10  BY 
WHICH  THE  GIVEN  VALUE  IS  TO  BE  MULTIPLIED. 
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TABLE  A- 3 

VALUES  OF  THE  DIFFUSION 
CONSTANT  FOR  18  GROUPS 


GROUP 

D (CM) 

D (CM) 

D (CM) 

NUMBER 

GRAPHITE 

FUEL 

GRAPHITE 

REGIONS  1 AND  5 

REGIONS  2 AND  4 

REGION  3 

1 

5.17803 

4.89794 

5.02004 

2 

4.25305 

3.92550 

3.92556 

3 

2.51475 

3.21046 

2.51475 

4 

1.88733 

2.66530 

1.88733 

5 

1.84187 

3.59117 

L. 14187 

6 

1.04608 

2.52128 

1.04045 

7 

0.97892 

2.08773 

0.97358 

8 

0.94844 

2.03009 

0.94303 

9 

0.96784 

1.56574 

0.96213 

1G 

0.95491 

1.67446 

0.94917 

11 

0.94745 

1.54715 

0.94185 

12 

0.94705 

1.76832 

0.94207 

13 

0.94772 

1.61025 

0.93559 

14 

0.93653 

1.45434 

0.93372 

15 

0.93525 

1.31442 

0.93242 

16 

0.93480 

1.19188 

0.93196 

17 

0.93481 

1.10335 

0.93197 

18 

0.9155 

0.3339 

0.9155 
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TABLE  A- 4 

VALUES  OF  THE  ABSORPTION 
CROSS-SECTION  FOR  18  GROUPS* 


GROUP  SIGMA  ABS.  (1/CM) 

NUMBER  FUEL 

REGIONS  2 AND  4 


1 

2 

3 

4 

5 

6 

7 

8 
9 

1C 

11 

12 

13 

14 

15 

16 

17 

18 


1 .96230(-3) 
7.177391-4) 
3.9C592(-4) 
2 . 52664 ( -4 ) 
2 . 12427 (-4 ) 
4.88448 (-4) 
3.7122K-4) 
1.056 10(-3) 
7 .492  5 1 ( -4 ) 
6.78011 (-4) 
6.82464 (-4 ) 
1. 130671-3 ) 
1. 502981-3) 
1.2C770(-3) 
3.41225{-3) 
5 . 57001 ( -3 ) 
9.15107(— 3) 
6.798  (-2) 


* THE  ABSORPTION  CROSS-SEC  I ION  FOR  GRAPHITE  IS  NEGLIGIBLY 
SMALL  FOR  ALL  GROUPS  EXCEPT  THE  18-TH.  FOR  GROUP  NUMBER  18 
THE  ABSORPTION  CROSS-SECTION  IS  0.0002414  (1/CM)  IN 
REGIONS  1,  3 AND  5. 


NOTE 


NUMBER  IN  PARENTHESES  INDICATES  THE  POWER  OF  10  BY 
WHICH  THE  GIVEN  VALUE  IS  TO  BE  MULTIPLIED. 
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TABLE  A- 5 

VALUES  OF  NU*S IGMA  FISSION 
AND  CHI  FOR  18  GROUPS 


GROUP  NC*  S IGMA  F (1/CM)  CHI 

NUMBER  FUEL 

REGIONS  2 AND  4 


1 

1.975051-3) 

0.02307 

2 

1.144571-3) 

0.10749 

3 

1.075981-3) 

0.20896 

4 

5.494401-4) 

0.52700 

5 

3.591601-4) 

0.10037 

6 

4.272901-4) 

0.01702 

7 

4.936401-4) 

0.00836 

8 

5.527301-4) 

0.00405 

9 

6.596601-4) 

0.00193 

10 

7.  721901-4) 

0.00092 

11 

9.201701-4) 

0.00044 

12 

1.031961-3) 

0.00021 

13 

1.215421-3) 

0.00010 

14 

1. 578191-3) 

0.C0007 

15 

2.937951-3) 

0.00002 

16 

5.821721-3) 

0.00000 

17 

1 . 14071 (-2) 

0.00000 

18 

1.108  (-1) 

0.00000 

NOTE  NUMBER  IN  PARENTHESES  INDICATES  THE  POKER  OF  10  BY 

WHICH  THE  GIVEN  VALUE  IS  TO  BE  MULTIPLIED. 
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TABLE  A-6 

VALUES  OF  THE  TOTAL  SCATTERING 
CROSS-SECTION  FOR  18  GROUPS 


GROUP 


FROM 

1 

TO 

I 

FROM 

1 

TO 

2 

FROM 

2 

TO 

2 

FROM 

1 

TO 

3 

FROM 

2 

Tu 

3 

FROM 

3 

TO 

3 

FROM 

1 

TO 

4 

FROM 

2 

TO 

4 

FROM 

3 

TO 

4 

FROM 

4 

TO 

4 

FROM 

1 

TO 

5 

FROM 

2 

TO 

5 

FROM 

3 

TO 

5 

FROM 

4 

TO 

5 

FROM 

5 

TO 

5 

FRUM 

1 

TO 

6 

FROM 

2 

TO 

6 

FROM 

3 

TU 

6 

FROM 

4 

TO 

6 

FROM 

5 

TO 

6 

FROM 

6 

TU 

6 

FROM 

1 

TO 

7 

FROM 

2 

TO 

7 

FROM 

3 

TO 

7 

FROM 

4 

TO 

7 

FROM 

5 

TO 

7 

FROM 

6 

TO 

7 

FROM 

7 

TO 

7 

FROM 

I 

TO 

8 

FROM 

2 

TO 

8 

FROM 

3 

TO 

8 

FROM 

4 

TU 

8 

NOTE  NUMBER 

WHICH 


SIGMA  S { 1/CM) 
GRAPHITE 
REG.  1 AND  5 

3.901531-2) 

3.158751-2) 

8 . 06 1 A 7 ( —2 ) 
7.683481-3) 
3.685851-2) 
8.721611-2) 

1.  181531-2) 

1 • 78686 ( -3 ) 
5.641911-2) 

1 . 806  76 ( - 1 ) 
2.090261-3) 
4.627311-4) 
O.OOOOO(-O) 
2.783301-2) 
2.563591-1 ) 
3.928171-4) 

1. 070591-4) 
O.OOOOO(-O) 
0.00000(-0 ) 
4.881481-2 ) 
2.333071-1 ) 
O.OOOOO(-O) 
O.OOOOO(-O) 
0.000001-0) 
O.OCOOO(-O) 
O.OOOOOI-O) 
1.080641-1 ) 
2.4593K-1  ) 
O.OOOOO(-O) 
O.OOOOO(-O) 
0.000001-0 ) 
O.OOOOO(-O) 


IN  PARENTHESES 
THE  GIVEN  VALUE 


SIGMA  S (1/CM) 
FUEL 

REG.  2 AND  4 

2.748991-2) 

1 .88118 (-2) 
4.208501-2) 
1.26735(-2) 
3.04768 (-2) 

6 . 53425 ( -2 ) 
2.288621 -2 ) 

3. 53009 (-2 ) 
6.99061 (-2) 
1.417491-1) 
4.622241-3) 
6.875501-3) 

1 .23346 (-2) 
5.63396 (-2) 

1 .91525 (-1 ) 
5.083581-4) 
1.137741-4) 
2.458041-3) 
1.266531-2) 
7.992511-2) 
2.022251-1) 
3.000151-4) 
6.796391-4) 
1.476021-3) 
7.418951-3) 
4.578861-2) 
1.467021-1) 
2.227811-1) 
1.788541-4) 
4.086711-4) 
8.907521-4) 
4.501111-3) 


INDICATES  THE 
IS  TO  BE  MULTI 


SIGMA  S 1 1/CM) 
GRAPHITE 
REGION  3 

4.386741-2) 
2.962911-2) 
7.917251-2 ) 
7.683481-3) 
3.685851-2) 
6.589291-2) 
1.181531-2  ) 
1.786861-3 ) 
5.641911-2) 
1.716161-1 ) 
2.145241-3  ) 
3.837321-4) 
C.00G001-0) 
2.783301-2) 
2.579121-1) 
3.868401-4) 
9.022811-5) 
0.000001-0) 

C. 000001-0 ) 
4.671881-2 ) 
2.347031-1) 
0.000001-0) 
0.000001-0) 
0.000001-0) 
O.OCOOO(-O) 
0.000001-0 ) 
1.066031-1 ) 
2.472631-1 ) 
C.OCOOO(-O) 
0.000001-0) 
0.000001-0) 
0.000001-0) 


POWER  OF  10  BY 
PLIED. 
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TABLE  A- 6 
VALUES  GF  THE  TOTAL 
CROSS-SECTION  FOR 


(CONTINUED) 
SCAT  TERING 
18  GROUPS 


GROUP 


SIGMA  S { 1 /CM ) SIGMA  S (1/CM)  SIGMA  S l 1/CM ) 


GRAPHI TE 

FUEL 

GRAPHITE 

REG.  1 AND  5 

REG.  2 AND  4 

REGION  3 

FROM 

5 

TO 

8 

o.ooooot 

-0) 

2.782431-2) 

o.ooooot-o 

FROM 

6 

TO 

8 

o.ooooot 

-0) 

8. 10891 (-2) 

O.OQOOOt-O 

FROM 

7 

TO 

8 

1. 15084( 

-1  ) 

1 .71677 (-1 ) 

1. 13640(-1 

FROM 

8 

TO 

8 

2 . 544  32 ( 

-1  ) 

1.7412K-1) 

2. 55954  t - 1 

FROM 

1 

TO 

9 

o.ooooot 

-0) 

1.072341-4) 

c.ooooot-o 

FRGM 

2 

TO 

9 

0.00000( 

-0) 

2 . 46426 ( -4 ) 

0.000001-0 

FROM 

3 

TO 

9 

o.ooooot 

-0) 

5 . 38374 ( -4 ) 

c.ooooot-o 

FRGM 

A 

TO 

9 

0.000001 

-0) 

2.73058(-3) 

o.ooooot-o 

FROM 

5 

TO 

9 

o.ooooot 

-0) 

1 .69001 t-2) 

o.ooooot-o 

FROM 

6 

TO 

9 

o.ooooot 

-0) 

4 . 92600 ( ~2 ) 

o.ooooot-o 

FROM 

7 

TO 

9 

o.ooooot 

-0) 

9.67055(-2) 

c.ooooot-o 

FROM 

8 

TO 

9 

1. 18731 ( 

-1 ) 

1.905471-1) 

1.172091-1 

FROM 

9 

TO 

9 

2.499401 

-1 ) 

2 . 7298  7 (-1 ) 

2.515001-1 

FROM 

2 

TO 

10 

o.ooooot 

-0) 

1 .49252 (-4 ) 

o.ooooot-o 

FROM 

3 

TO 

10 

o.ooooot 

-0) 

3 . 26  556  t -4 ) 

o.ooooot-o 

FROM 

4 

TO 

10 

o.ooooot 

-0) 

1 . 66030 ( -3 ) 

c.ooooot-o 

FROM 

5 

TO 

10 

o.ooooot 

-0) 

1 .02846 (-2 ) 

o.ooooot-o 

FROM 

6 

TO 

10 

o.ooooot 

-0) 

2 . 99806 (-2 ) 

c.ooooot-o 

FROM 

7 

TO 

10 

o.ooooot 

-0  ) 

5 . 88598 ( -2 ) 

o.ooooot-o 

FROM 

8 

TO 

10 

o.ooooot 

-0) 

1. 102 1 0 ( -1 ) 

c.ooooot-o 

FROM 

9 

TO 

10 

1.16811  ( 

-1 ) 

2.0871K-1) 

1. 152331-1 

FROM 

10 

TO 

10 

2 . 526  7 1 ( 

-1) 

1 .67039(-l) 

2 . 54247 ( -1 

FROM 

3 

TO 

11 

o.ooooot 

-0) 

1.955011-4) 

0.00000(-0 

FROM 

4 

TO 

1 1 

o.ooooot 

-0) 

9 . 9 548  3 ( -4 ) 

c.ooooot-o 

FROM 

5 

TO 

11 

o.ooooot 

-0) 

6.169711-3) 

o.ooooot-o 

FROM 

6 

TO 

11 

o.ooooot 

-0) 

1. 798661-2) 

o.ocooot-o 

FROM 

7 

TO 

11 

o.ooooot 

-0  ) 

3 . 53 1 36 ( -2 ) 

C.0G0001-0 

FROM 

8 

TO 

1 1 

o.ooooot 

-0) 

6.612331-2) 

c.ooooot-o 

FROM 

9 

TO 

11 

o.ooooot 

-0) 

1 . 2 3482  t - 1 ) 

o.ocooot-o 

FROM 

10 

TO 

11 

1 . 1 809  1 ( 

-1) 

2 . 18606 (-1 ) 

1. 16503(-1 

FROM 

11 

TO 

1 1 

2 . 5474  7 ( 

-1 ) 

1.866861-1 ) 

2.563161-1 

FROM 

4 

TO 

12 

o.ooooot 

-0  ) 

6 . 07499 ( -4 ) 

C.00000(-0 

NOTE 

NUMBER 

IN  PARENTHESES 

INDICATES  THE 

POWER  OF  10 

WHICH 

THE  GIVEN 

VALUE 

IS  TO  BE  MULTIPLIED. 
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TABLE  A-6  ( CON  T I IMUED  ) 
VALUES  UF  THE  TOTAL  SCATTERING 
CRUSS-SECTI ON  FOR  18  GROUPS 


GROUP  SIGMA  S (1/CM)  SIGMA  S ( 1/CM ) SIGMA  S (1/CM) 


GRAPHITE 

FUEL 

GRAPHI TE 

REG.  1 AND  5 

REG.  2 AND  4 

REGION  3 

FROM 

5 

TO 

12 

O.COQCO(-O) 

3.766431-3) 

O.OCOOO(-O) 

FROM 

6 

TO 

12 

O.OOOOO(-O) 

1.09808 (-2) 

C.COOOO(-O) 

FROM 

7 

TO 

12 

O.OOOOO(-G) 

2 . 15594 ( -2 ) 

O.OOOOO(-O) 

FROM 

8 

TO 

12 

O.OOOOO(-O) 

4.03698 ( -2 ) 

O.OOOOO(-O) 

FROM 

9 

TO 

12 

O.OOOCO(-O) 

7 . 53899 ( -2 ) 

O.OCOOO(-O) 

FROM 

10 

TO 

12 

O.OOOOO(-O) 

1.3 1674 (-1 ) 

C.OCOOO(-O) 

FROM 

11 

Tu 

12 

1.  190231-1  ) 

2.302151-1 ) 

1. 17  444 ( - 1 ) 

FROM 

12 

TO 

12 

2 * 56700 ( - 1 ) 

1 . 9432  8 ( - 1 ) 

2. 5752 1 ( -1 ) 

FROM 

5 

TO 

13 

O.OOOCO(-O) 

2.276471-3) 

C.OCOOO(-O) 

FROM 

6 

TO 

13 

O.OOOOO(-O) 

6.63709 (-3) 

C.OOOOO(-O) 

FROM 

7 

TO 

13 

O.OOOOO(-O) 

1 . 30  3 1 3 ( -2  ) 

C.OOOOO(-O) 

FROM 

8 

TO 

13 

O.OOOOO(-O) 

2.440131-2) 

O.OCOOO(-O) 

FROM 

9 

TO 

13 

O.OOOCO(-O) 

4 . 55693 ( -2 ) 

C.OCOOO(-O) 

FROM 

10 

TO 

13 

0.C00C0(-0 ) 

7.959091-2) 

C.OCOOO(-O) 

FROM 

11 

TO 

13 

O.OOOOO(-G) 

1 . 36850 ( -1 ) 

O.OOOOO(-O) 

FROM 

12 

TO 

13 

1.  1879  1 { - 1 ) 

2 . 40694 ( — 1 ) 

1.1 7966 ( - 1 ) 

FROM 

13 

TO 

13 

2.57314(-1 ) 

i.99820(-l) 

2.581241-1 ) 

FROM 

6 

TU 

14 

O.COOOO(-O) 

6 . 5 1449 ( -3 ) 

C.OOOOO(-O) 

FROM 

7 

TO 

14 

O.OOOOO(-O) 

1 . 27907 ( -2 ) 

O.OOOOO(-O) 

FRUM 

8 

TO 

14 

O.OOOOO(-O) 

2 • 395 10 ( -2 ) 

O.OCOOO(-O) 

FROM 

9 

TO 

14 

O.OOOOO(-O) 

4.47287 (-2 ) 

O.OOOOO(-O) 

FROM 

10 

TU 

14 

O.COOOO(-O) 

7 . 8 12  36 ( -2 ) 

C.OOOOO(-O) 

FROM 

11 

TO 

14 

O.OOOOO(-O) 

1 . 34328 ( -1  ) 

0. 00000 ( -0 ) 

FROM 

12 

TO 

14 

0.00000(-0 ) 

2. 31363 (-1 ) 

O.OOOOO(-O) 

FROM 

13 

TO 

14 

1 . 19Q59(- 1 ) 

3.90856{-l) 

1.18247 (-1 ) 

FROM 

14 

TO 

14 

3 . 172 1 1 ( — 1 ) 

3.21039  (-1 ) 

3. 18534 ( -1 ) 

FROM 

7 

TO 

15 

O.OOOOO(-O) 

4.659791-3) 

C.00000(-0 ) 

FROM 

8 

TO 

15 

O.OOOOO(-O) 

8. 72568(-3) 

O.OCOOO(-O) 

FROM 

9 

TO 

15 

O.OOOCO(-O) 

1 . 6295  3 ( -2 ) 

O.OCOOO(-O) 

FROM 

10 

TO 

15 

O.COOOO(-0 ) 

2 • 846 18 (-2 ) 

O.OCOOO(-O) 

NOTE 

NUMBER 

IN  PARENTHESES 

INDICATES  THE 

POWER  OF  10  BY 

WHICH  THE  GIVEN  VALUE  IS  TO  BE  MULTIPLIED. 
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TABLE  A-6  { CON  T I NUED ) 
VALUES  OF  THE  TOTAL  SCATTERING 
CROSS-SECTION  FOR  18  GROUPS 


GROUP 


FROM 

11 

TO 

15 

FROM 

12 

TO 

15 

FROM 

13 

TO 

15 

FROM 

14 

TO 

15 

FROM 

15 

TO 

15 

FROM 

8 

TO 

16 

FROM 

9 

TO 

16 

FRUM 

10 

TO 

16 

FROM 

11 

TO 

16 

FROM 

12 

TO 

16 

FROM 

13 

TO 

16 

FROM 

14 

TO 

16 

FROM 

15 

TO 

16 

FROM 

16 

TO 

16 

FROM 

9 

TO 

17 

FROM 

10 

TO 

17 

FROM 

11 

TO 

17 

FROM 

12 

TO 

17 

FRUM 

13 

Tu 

17 

FRUM 

14 

TO 

17 

FROM 

15 

TO 

17 

FROM 

16 

TO 

17 

FROM 

17 

TO 

17 

FROM 

10 

TO 

18 

FROM 

11 

TO 

18 

FROM 

12 

TO 

18 

FROM 

13 

TO 

18 

FRUM 

14 

TO 

18 

FROM 

15 

TO 

18 

FROM 

16 

TO 

18 

FRUM 

17 

TO 

18 

SIGMA  S (1/CM) 
GRAPHITE 
REG.  1 AND  5 

O.OOOOO(-O) 
O.OOOOO(-U) 
O.OOOCO(-O) 
6.0C05  1 (-2 ) 
3.176861-1 ) 
O.OOOOO(-O) 
O.OCOOO(-0 ) 
O.OOOOO(-O) 
O.OOOOO(-O) 
O.OOOOO(-G) 
0.000001-0 ) 
O.OOOOO(-O) 

6. 009221-2) 

3 . 47573 ( - 1 ) 
O.OOOOO(-O) 
O.OOOOO(-O) 
O.OOOOO(-O) 
O.OOOOO(-O) 
O.OOOOO(-O) 
O.OOOOO(-O) 
O.OOOOO(-O) 

3 . 04048 ( -2 ) 

3. 640  78 ( - 1 ) 
O.OOOOO(-O) 
O.OOOOO(-O) 
O.OOOOO(-O) 
C.OOOOO(-O) 
O.OOOOO(-O) 
O.OOOOO(-O) 
O.OOOOO(-O) 
1.38991 (-2 ) 


SIGMA  S ( 1/CM) 
FUEL 

REG.  2 AND  4 

4.89381 (-2) 
8.42902 (-2) 

1 . 40782 ( - 1 ) 
3.150451-1) 
3.241901-1 ) 
4.41453 (-3) 

8 . 24424 ( -3 ) 

1 .43996 (-2 ) 
2.47591 (-2) 

4 • 26447 ( -2 ) 

7. 12258(— 2) 

1 . 5 8284 {- 1 ) 

4 . 28848 (-1 ) 
4.892561-1) 
1.27500 (-3) 

2 . 22695 (-3 ) 

3 . 82909 ( -0 ) 

6 . 595 1 6 ( -3 ) 

1. 10153 (—2 1 
2 .44792 (-2 ) 

6 . 59962 ( -2 ) 
3.34482 (-1 ) 
6.62588(-l) 

2 .50141 (-5) 
4.3010Q(-5) 
7.407971-5) 
1.237291-4) 

2 . 7496  1 (-4 ) 
7.412991-4) 
3.7  4469 { -3 ) 

1 • 71 589 (—1 ) 


SIGMA  S (1/CM) 
GRAPHITE 
REGION  3 

O.OOOOO(-O) 
O.OOOOO(-O) 
C.OOOOO(-O) 
5.867781 —2 ) 

3. 19007 (-1 ) 
C.OOOOO(-O) 
O.OOOOO(-O) 
O.OCOOO(-O) 
0.000001-0) 
O.OOOOO(-O) 
O.OOOOO(-O) 
O.OCOOO(-O) 

5 . 87692 ( -2 ) 
3.491491-1 ) 
O.OOOOO(-O) 
O.OCOOO(-O) 
O.OOOOO(-O) 
C.OOOOO(-O) 
C.OOOOO(-O) 
0.00000(-0  ) 
O.OOOOO(-O) 

2 . 88290 ( -2 ) 

3. 65784 (-1) 
O.OOOOO(-O) 
O.OOOOO(-O) 
O.OOOOO(-O) 
0.00000(-0 ) 
C.OOOOO(-O) 
O.OOOOO(-O) 
O.OOOOO(-O) 

1 . 2 1929 ( -2 ) 


NOTE  NUMBER  IN  PARENTHESES  INDICATES  THE  POWER  OF  10  BY 

WHICH  THE  GIVEN  VALUE  IS  TO  BE  MULTIPLIED. 
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TABLE  A- 7 

AVERAGE  INVERSE  SPEED 
FOR  18  GROUPS 


GROUP 


L/V  (SEC/CM) 


1/V  ( SEC /CP ) 


1/V  (SEC/CP) 


NUMBER 

GRAPHITE 

FUEL 

GRAPHITE 

REG.  1 AND 

5 

REG.  2 AND 

4 

REGION  3 

I 

2.765161- 

10) 

2.692351- 

10) 

2.721531-10) 

2 

3.342391- 

10) 

3.414061- 

10) 

3.38431 (-10) 

3 

4.384381- 

10) 

4.338191- 

10) 

4 . 39308 ( - 10 ) 

4 

7.638971- 

10) 

7.305221- 

10) 

7. 59417 (-10) 

5 

1.354261- 

9) 

1.331991- 

9) 

1.348201-  9) 

6 

1.933971- 

9) 

1.932491- 

9) 

1.932131-  9) 

7 

2.479991- 

9) 

2.488891- 

9) 

2 .47783 ( - 9) 

8 

3.188531- 

9) 

3.177771- 

9) 

3.185641-  9) 

9 

4.093361- 

9) 

4.109761- 

9) 

4.089451-  9) 

LO 

5.25631 (- 

9) 

5.246891- 

9) 

5.251361-  9) 

II 

6. 74941 ( - 

9) 

6.739591- 

9) 

6.74  308 ( - 9) 

12 

8.662811- 

9) 

8.657651- 

9) 

8.658591-  9) 

13 

1.112351- 

8) 

1.111861- 

8) 

1.11 182 (—  8) 

14 

1.657391- 

8) 

1.655551- 

8) 

1 .65352 (-  8) 

15 

2.73271 (— 

8) 

2.729371- 

8) 

2 .72633 ( - 8) 

16 

6 . 3406 1 ( - 

8 ) 

6.307961- 

8) 

6 .28346 ( - 8) 

17 

4.902121- 

7) 

4.777791- 

7) 

4 . 74673 ( - 7) 

18 

5. 053511- 

6) 

4.907651- 

6) 

5.05351 (-  6) 

NOTE 

NUMBER  IN 

PARENTHESES  INDICATES 

THE  POWER  OF  10  B 

WHICH  THE  GIVEN  VALUE  IS  TO  BE  MULTIPLIED. 
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TABLE  A- 8 
AVERAGE  ENERGY 
FOR  18  GROUPS 


GROUP 

NUMBER 

E (MEV) 
GRAPHITE 
REG.  1 AND  5 

E (MEV) 
FUEL 

REG.  2 AND 

4 

E (MEV) 
GRAPHITE 
REGION  3 

1 

6 .89044 ( - 0) 

7 .32742 (— 

0) 

7 . 15225  ( - 

2 

4 * 76849 ( — 0) 

4 . 56614 ( - 

0 ) 

4 . 65050 ( - 

3 

2.763651-  0) 

2.83051 (- 

0) 

2.75700(- 

A 

1.025781-  0) 

1. 12520(- 

0) 

1 .03842 (— 

5 

3.043521-  1) 

3 . 14496  ( - 

1 ) 

3 . 07 l 30 ( - 

6 

1 .42451 ( - 1) 

1.42671 (- 

1 ) 

1 .42724 ( - 

7 

8.663081-  2) 

8.60053 (- 

2 ) 

8.678271- 

8 

5 . 24G64 ( - 2) 

5 . 2 7636 ( — 

2) 

5.25023 (— 

9 

3.179841-  2) 

3 . 1 54  14  ( - 

2 ) 

3. 18598 (- 

10 

1.928431-  2) 

1.935401- 

2) 

1 .93209 ( - 

1 1 

1 . 16959 { - 2) 

1. 17303(- 

2 ) 

1 . 1 7 1 81 ( — 

L2 

7.099881-  3) 

7 . 10841 {— 

3) 

7.1 0686  ( - 

13 

4.306121-  3) 

4 .30993 ( - 

3) 

4.310251- 

1 A 

2.031421-  3) 

2.03609(- 

3) 

2.04120(- 

15 

7.472421-  4) 

7.49105{— 

4) 

7 . 50832  ( - 

16 

1.657931-  4) 

1.67736(- 

4) 

l .69125 (— 

17 

6.250361-  6) 

6 . 7 1803 ( - 

6) 

6.879221- 

18 

4.737351-  8) 

5 . 502  30 ( - 

8) 

4 . 7 37  35 ( - 

NOTE 

NUMBER  IN  PARENTHESES  INDICATES 
WHICH  THE  GIVEN  VALUE  IS  TO  BE 

THE  POWER  OF  10 
MULTIPLIED. 
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